大規模海洋數據同化的并行優化

蔡迪 洪學海 肖俊敏 譚光明

（中國科學院計算技術研究所 北京 100190）

海洋數據同化是一種從氣候數值預測領域發展而來的，將物理模型與觀測資料相結合的預測校正方法[1-3]，廣泛運用于大氣[4]、海洋[5]和地表[6]等諸多領域.目前最為常見的數據同化方法分別為集合卡爾曼濾波算法（FnKF）[7]、集合最優插值算法[8]，以及三維、四維變分方法等[9-12].其中集合最優插值由集合卡爾曼濾波發展而來.

隨著航海事業的不斷發展，捕魚、軍事、資源開采和生物研究等海洋相關的活動日益增多.這些日益增多的海洋活動對海洋預報的準確性和實時性提出了更高的要求，因此許多準實時的海洋觀測系統被開發出來[13].隨著海洋模式的不斷發展，海洋數據同化的分辨率也不斷提高.在其他因素不變的情況下，水平方向的分辨率每提高10 倍，對應的計算量和數據量將提高100 倍，這大大增加了同化系統對內存、I/O 以及計算能力的要求.海量的數據讀取和計算任務對系統的實時性帶來了很大的挑戰.我國自主研發的天河2 號超算平臺憑借其優越的I/O 性能以及計算的性能曾多次登頂世界超算的榜首，可以很大程度上滿足同化程序對于內存、I/O 以及計算性能的需求，因此許多大數據學者都以天河2 號作為研究和開發平臺[14-24].如何開發一個高效的并行算法提取出海洋數據同化的特性，并與天河的存儲以及執行架構相結合，最大程度地提升天河系統的內存、I/O 以及計算資源的利用率顯得尤為重要.

本文算法主要針對基于中國科學院大氣物理研究所（Institute of Atmospheric Physics，Chinese Academy of Sciences，IAP）、大氣科學和地球流體力學數值模擬國家重點實驗室（State Key Laboratory Modelling for Atmospheric Sciences and Geophysical Fluid Dynamics，LASG）發展的 LASG/IAP 氣候系 統海洋模式（LASG/IAP climate ocean model，LICOM）[25-27]，在局部集合最優插值的算法框架下，對大規模高分辨率的海洋數據同化程序進行并行優化.由于本文采用的是高分辨率的海洋觀測數據，因此同化時存在計算量大、內存需求高、I/O 時間長等問題，在文獻[28]中已經針對此這個模式下的數據同化程序做了一些優化工作，取得了不錯的效果，但是其優化時沒有深入理解海洋數據同化的數據特性、計算特性，沒有考慮所使用的超算平臺架構，以及沒有考慮算法的時間局部性和空間局部性的問題，因此仍存在很大的優化空間.本文在文獻[28]的基礎上，用其優化后的程序作為本文的對比程序（下文稱之為I/O 代理程序），并進一步將海洋數據同化的數據特性、計算特性與所使用的超算平臺的架構特性相結合，并結合時間局部性和空間局部性，設計了一款高效的大規模海洋數據同化并行算法.

1 算法及程序流程

1.1 同化方法

局部集合最優插值算法是由集合最優插值進行可并行性優化后發展而來[29]，具有更好的并行性.

本文局部集合最優插值的實現過程表示為

其中，φα為分析場數據，φ為背景場數據，α為一個在區間（0，1]上的隨機因子；A為靜態樣本數據，A′為擾動預測矩陣，由A中心化而來，A′A′T為背景誤差協方差的估計；H為測量算子，代表真實模型與預測模型之間的映射；γ為擾動測量誤差向量，γγT為觀測誤差協方差；d′=d-Hφ為觀測增益，其中d為擾動測量向量.

令

由于集合擾動誤差和測量誤差無關，所以由式（2）得

結合式（1）可得

SVD 分解以及式（7）便是本文各個海洋格點的主要計算部分，其中最為耗時的便是SVD 分解過程.在SVD 分解完成之后，將式（7）分為5 個步驟計算：

在局部最優插值的算法中，整個同化數據區域按照經緯度被劃分成很多個格點，對于每一個格點只要可以獲取到周圍（2×prep_rx）×（2×prep_ry）范圍內的數據便可以對該節點進行數據同化，其中在本文的分辨率下prep_rx=105，prep_ry=19. 這樣良好的并行特性十分便于計算的并行化，但是這也導致了每個格點在數據分發時需要帶邊數據，在進行通信時會產生大量的通信量，增加數據通信的時間，大大降低程序的并行效率. 因此為了提升程序的并行效率、提高并行系統的帶寬利用率以及減少通信時間，提升程序的并行效率尤為重要，這些將在下文做詳細介紹.

1.2 程序流程

文獻[28]的程序的流程圖如圖1 所示. 并本文以此作為對比程序，并在此基礎上進行優化.

1.3 存在的挑戰

1）負載均衡的挑戰

Fig.1 Flowchart of LICOM data assimilation pragram based on the I/O agents strategy圖1 基于I/O 代理策略的LICOM 數據同化程序流程圖

由圖2 可以看出，I/O 代理程序各個進程的計算負載極不均衡，大量的進程沒有計算任務，大部分計算任務被分配到了少量的進程上，且這些進程內部負載也十分的不均衡，存在極大的進程浪費. 在I/O代理程序中，同化程序以同化數據中海洋格點的數目為負載均衡的依據，對進程進行計算任務的劃分.這樣的劃分方式只是從同化流程的角度理解海洋數據同化，沒有深入理解同化算法以及進程的計算過程. 因此在這種負載均衡方式下，雖然各個進程表面上看起來負載均衡，實際上執行時存在很大的負載不均衡性. 雖然每個進程獲得的可同化格點的數目是一樣的，但并不是每個海洋格點都是可以同化的，當格點周圍的nobs≤2 時（其中nobs代表同化格點周圍有效的觀測數據點的數量），此部分格點是不可以同化的. 最壞的情況是，一些進程獲取的海洋格點全部需要同化，而一些海洋格點獲取的海洋格點完全不需要同化， 這要導致極大的負載不均衡性，從而造成大量計算節點的浪費.

2）讀取的挑戰

由圖3 可以看出I/O 代理程序在讀取64.8 GB 的數據時較為耗時，并行讀取效率較低. I/O 代理程序讀取時采用代理節點對數據進行讀取，之后由代理節點對計算進程進行數據分發. 這樣的設計方式僅在上層算法層面對同化程序的讀取進行了優化，并未考慮到超算集群的文件系統以及通信架構，因此并未充分挖掘Lustre 并行文件系統的并行讀取能力，并行帶寬利用率較低.

Fig.2 The assimilation time of processes under 1 000 cores for the I/O agents program圖 2 I/O 代理1 000 核下各個進程的同化時間

Fig.3 The time for the I/O agents program to read 64.8 GB of data under different numbers of cores圖 3 不同核數下I/O 代理程序讀取64.8 GB 數據的時間

3）通信的挑戰

由圖4 可以看出I/O 代理程序的通信不僅耗時較長且隨著核數增加顯著增長.正如讀取的挑戰所述，I/O 代理程序并未考慮超算集群的節點通信架構，只是簡單地調用函數MPI_Scatter進行數據分發，存在很大的通信負載不均衡，超算集群通信通道利用率非常低.

Fig.4 Communication time of I/O agents program under different numbers of cores圖 4 不同核數下I/O 代理程序的通信時間

4）寫回的挑戰

由圖5 可以看出I/O 代理程序的寫回十分的耗時，且隨著核數增加十分的不穩定.I/O 代理程序寫回時采用各個進程直接寫回對應數據塊的策略，這樣的策略會導致海量的進程同時密集地訪問同一個對象存儲目標（object storage target，OST），而OST 并不能同一時間處理這么多I/O 進程的請求，造成極大的沖突，特別是核數達到千以上時.此外海量的進程同時訪問同一塊磁盤也會導致嚴重的磁頭爭用和海量的尋址.

Fig.5 Write back time of I/O agents program under different cores圖 5 不同核數下I/O 代理程序的寫回時間

1.4 本文貢獻

本文針對1.3 節所述的I/O 代理程序存在的4 個挑戰，通過分析同化程序整體的算法流程，對算法進行重組并結合超算集群的存儲架構、節點分布特性以及通信特性，提出了基于計算拓撲圖的2 層負載均衡策略以及基于Lustre 文件存儲架構和超算集群特性的并行優化策略，并在此基礎上進一步提出了計算、讀取通信、寫回3 層重疊策略.本文的總技術路線如圖6 所示.

2 負載均衡優化

2.1 基于計算拓撲圖的負載均衡優化

Fig.6 General technical route圖 6 總體技術路線

本文對海洋數據同化算法和進程實際的計算流程進行深入的分析，觀察到nobs才是各個進程實現負載均衡的關鍵.而且實際同化過程中最耗時的部分是SVD 分解過程，其計算復雜度為nobs3，其中nobs的取值范圍為3～1 000.由此也可以看出基于海洋格點的負載均衡方法中，有的海洋格點nobs很大，有的海洋格點nobs很小，計算復雜度最大相差1 0003/27，為370 多萬倍，這樣的負載均衡方法十分不合理.本文算法針對實際執行流程中的這一特性，在實際執行同化前加入了一個預同化步驟，分別計算出每一個海洋格點進行數據同化的計算復雜度，在原本海洋格點分布圖的基礎上生成基于計算復雜度的計算拓撲圖，并以此作為負載均衡劃分的依據，在最大程度上實現了各個進程的負載均衡.不僅如此，實現的計算拓撲圖還可以很大程度上反映出本次數據同化計算任務的分布特性，顯示出整個海洋同化數據中的對應海洋格點是否需要同化，以此實現了針對海洋同化數據計算分布特性的讀取寫回以及通信優化，這部分將在2.2 節進行介紹.

2.2 2 層負載均衡

Fig.7 Two-layer load balancing圖 7 2 層負載均衡

如圖7 所示，假設整個海洋同化數據的計算量為90，每個節點（Split 組）內5 個核，其中0 號核為該Split 組的Master.生成計算拓撲圖后，為了方便展示，假設將整個同化數據分為4 組（25，20，30，15）并計算出每個組的計算量.本文先在組的維度按照經度進行1 層的負載均衡，以第1 組為例，其計算量為25，本文算法為這個組分配5 個Master，由于地理區域上的觀測點分布是不均勻的，因此每個Master 負責的地理區域大小可能有大有小，每個Master 負責的計算量為5.在每個Split 組內部，本文算法再進行了1層的負載均衡，將計算量再次平均劃分給節點內的計算進程.這里的Split 組屬于不同的節點，分別使用節點間和節點內通信通道，這樣2 層的負載均衡方式不僅可以最大限度地實現計算進程的負載均衡，而且更為第3 節的讀取和通信分組做了鋪墊，第3 節的讀取和通信優化策略也充分考慮節點間和節點內的通信通道問題，并根據負載均衡的分組，實現了相應的讀取和通信的優化策略.通過這樣的分組方式，每個組的Reader 對Master 的數據分發，以及Master對Split 組的數據分發分別使用的是節點間和節點內的通信通道，互不影響，可以極大地提升通信效率.

3 基于Lustre 文件存儲架構和超算集群特性的并行優化

3.1 基于OST 的讀取優化

基于LICOM 模式的數據同化本身涉及大量的文件讀取，在同化時間沒有充分優化時，讀取時間可以忽略不計，例如在文獻[28]的Fortran 版本中，程序千核的運行時間為8 000 s 左右，而讀取時間占總時間的比重很小，可以忽略不計.但在同化時間充分優化后，讀取時間將在程序整體運行時間中占據很大的比重，特別是在本文負載均衡優化后，千核下程序其他部分的執行時間只有幾十秒，而程序讀取64.8 GB的數據需要10 多秒，因此讀取時間的優化也是提升并行同化程序性能的關鍵.

在傳統的讀取算法優化中，大多數的優化算法都是通過將數據重組，保證數據的成塊連續讀取，減少文件的尋址次數和I/O 的請求次數，從而減少文件的讀取時間.MPI-IO 是一種讀取通信優化策略，采用的就是這個思想，許多讀取密集型的算法以此作為讀取的優化策略，并在MPI-IO 的基礎上進行了優化[30-34].

傳統的讀取優化大多是在算法層面對I/O 進行優化，并沒有考慮過讀取進程是如何和底層文件存儲系統進行交互的，也沒有考慮超算集群的特性.因此這些優化算法往往只能一定程度上優化程序的讀取時間，如果舊的讀取策略已經是成塊連續讀取的，新優化算法的性能提升將十分的有限，而且隨著平臺網絡環境的不同有很大的波動.本文提出的讀取優化，將超算平臺的文件存儲架構和超算集群的特性結合起來，充分挖掘了Lustre 文件存儲架構的I/O帶寬和超算集群的節點通信通道的分配特點[35-42]，程序讀取64.8 GB 文件的讀取時間縮小到1.6 s 左右，I/O 性能顯著提升.

本文程序所使用的天河2 號超算平臺使用的文件系統為Lustre 文件系統.文件在存儲時系統按照分片的方式，將文件存儲到不同的OST 中，其中對象存儲服務器（object storage server，OSS）負責處理上層的I/O 讀寫請求以及調度OST.在忽略OSS 調度時間的影響下，如果OST 的數目為m，且OST 的I/O 帶寬為n（單位為GBps），那么理論上整個超算平臺的帶寬可以達到mn（單位為GBps）.在天河2 號超算平臺上，OST 的數目為12，OSS 的數目為2.如果忽略掉OSS的調度時間，天河系統的整體帶寬可達12n（單位為GBps）.在實測時，I/O 代理程序的吞吐量較低，遠遠沒有達到天河的峰值吞吐量.

首先，如圖8，9 所示，I/O 代理程序的I/O 讀取算法存在大量的讀取進程同時訪問少量OST 的情況，并且由于未對文件存儲方式進行設置，默認情況下同一個文件可能無序地存在幾個不同的OST，這就導致了雖然讀取的是不同的文件但不同的進程可能會訪問相同的OST.這樣的讀取方式使得同一時間

只有少量OST 超負荷運行，其他OST 則一直空閑，沒有充分利用OST，而且同時對少量OST 進行密集訪問會存在大量的I/O 請求沖突，增加了OSS 的調度難度，大大增加了OSS 的調度時間.在本文的算法中，對于64.8 GB 也就是12 個文件，每個文件5.4 GB 的數據，本文通過Lustre 的用戶接口，將每個文件分別存到不同的OST 中.在文件讀取時，本文只指定12個進程作為Reader 進行文件讀取.對于這12 個進程本文也做了特別的設定，本文將這12 個進程分別綁定到不同的天河節點上.2.2 節中提到，天河上的節點在通信時分別使用節點間和節點內的通信通道.

Fig.8 Reading strategy for the I/O agents pragram圖 8 I/O 代理程序的讀取策略

Fig.9 Reading strategy based on OST圖 9 基于OST 的讀取策略

如圖10 所示，天河上每個節點有24 個核，每個節點獨自占用1 個節點間通信通道，節點內部則走節點內的通信通道.這樣的設計方式可以充分利用每個節點的通信通道的帶寬，最大限度地避免對OST 請求的沖突.當OSS 數目足夠多時，每個讀取進程和OST 的I/O 請求完全并行，不受其他讀取進程的影響，實現完全并行的通信.但是正如前文提到的天河的OSS 只有2 個，因此無法實現理論上I/O 請求完全沒有沖突的情況，但是這種設計方式也可以最大限度地避免I/O 請求的沖突，并且同一時間充分利用天河可用的OST，實現完全并行讀取.

Fig.10 The layout and communication mode of supercomputing cluster computing node圖 10 超算集群計算節點布局與通信方式

3.2 3 層通信

1）3 層通信算法

本文的算法在讀取時設置了12 個Reader 對12個文件進行完全并行的讀取，接下來要解決Reader對讀取到的數據進行分發的問題.在I/O 代理程序中，讀取進程讀取到相應的數據之后，就直接對所有組內的進程進行數據分發.這樣的實現方式雖然簡單，但是效率卻十分的低下.在實際執行時1 個讀取進程可能要同時將數據分發給海量的進程，雖然MPI 在底層優化時，采用二分的方法進行數據分發的優化，但是并沒有關注圖10 中超算集群上節點內和節點間的通信通道問題.這樣可能存在一種情況，即在MPI_Scatter優化算法中，分發進程集中到1 個或幾個節點上，同時使用1 個節點的通信通道向外分發數據，這樣就會造成很大的網絡擁塞，分發的效率十分的低下.MPI_Scatter只是一種通用性的方法，實際執行的效果不佳.

Fig.11 Three layers of communication圖 11 3 層通信

本文根據圖10 所示的超算集群上節點間和節點內通信通道的配置特點設計了3 層通信算法.如圖11所示，首先設置12 個進程作為Reader 與OST 進行通信，獲取到相應的數據，其中每個Reader 分別分配到不同的節點上，獨立使用自己的通信通道與自己對應的OST 進行交互，與其他Reader 互不影響實現并行的OST 交互.其次在Reader 獲取相應的數據之后，將計算進程按照節點為單位分為Split 組，其中每個節點的0 號進程為這個Spilt 組的Master.0 號進程負責節點間的收發數據，獨占當前節點的通信通道.最后為了并行處理，本文將同化數據按照緯度ny方向橫向切塊，分成了group個Master 組.每個Master 組都需要12 個文件中對應的1 塊數據塊，因此本文為每個Master 組配置了12 個Reader，每輪由12 個Reader 去并行讀取12 個文件的對應數據塊，并行分發給對應的Master.這樣的設計方式可以將Reader和Master 間的通信集中到節點間的通信通道上，充分利用各個節點通信通道的同時，實現和下文節點內通信的隔離.在Master 獲取到相應的數據后，各個Master 通過節點內的通信通道向自己負責的Split 組分發數據.這樣設計方式可以實現節點間和節點內在不同通信通道的并行通信.

2）基于環的分發策略

每個Master 組的12 個Reader 對各自的Master 進行分發數據時，如果直接采用數據分發的方式，同一時間12 個Reader 同時對所有Master 發起通信，可能會出現“hot spot”，導致嚴重的通信擁塞.因此本文設計了圖12 所示的基于環的分發策略.

Fig.12 Ring-based distribution strategy圖 12 基于環的分發策略

本文首先將每個Master 組內的Master 再次平均分組，生成Scatter_Master 組，初始情況下每個Reader對應1 個Scatter_Master 組，負責向這部分的Master分發數據.在1 輪分發完成后，每個Reader 進行輪轉，如圖12 所示，依次循環直到所有的數據分發完畢.這樣的設計方式可以實現Reader 間完全并行的數據分發，不會產生“hot spot”和Reader 間的通信擁塞，極大地提升了通信的效率.

3.3 預同化優化和寫回優化

1）預同化優化

第2 節提到的基于計算拓撲圖的負載均衡策略中需要進行一步預同化的操作，這部分需要先讀取背景場的數據，然后根據背景場、觀測場和網格數據進行預同化，計算出每個海洋格點的計算復雜度，以此生成計算拓撲圖.這個過程涉及一張背景場文件讀取，這個背景場文件的大小為5.4 GB.在1 000 核下，這部分的讀取和通信時間為6 s 左右，且隨著核數呈線性增長.本文算法在4 000 核下優化到最后的整體運行時間為34 s 左右，因此預同化部分的優化也十分的重要.

本文將5.4 GB 的文件分別等分到group（group＜12）個OST 中，然后每個Master 組設置1 個Reader 對自己負責的OST 進行并行讀取.在讀取完成后仿照3.2 節提到的3 層通信，將數據分發給不同的組，各個組收到數據后進行預同化過程.

2）寫回優化

I/O 代理程序中，寫回時間占I/O 代理程序運行時間的27%，其時間為200～300 s 之間，且隨著超算集群用戶數量和網絡擁塞情況出現很大的波動.有時寫回時間甚至大于1500 s，不僅十分的耗時且十分的不穩定.這主要是由于I/O 代理程序采用各個計算進程獨自寫回所同化的數據，這樣的寫回策略會對OSS 發起數百萬次的寫回請求，同一時間OSS 需要處理幾千次的文件寫回請求，這對OSS 的調度和超算集群的網絡通信產生非常大的壓力，當超算集群用戶增多、通信網絡變得擁堵時，寫回的性能會受很大的影響.本文根據Lustre 文件存儲架構的特性，將寫回的文件分別根據Master 組等分到不同的OST 中，然后根據3 層通信，實現了3 層的gather 操作，將每個計算進程計算完成的數據經由Split→Master→Reader 收集到每個Master 組的Reader 0 中，由每個組的Reader 0 負責集中的寫回，將對應的數據寫回對應的OST 中，實現完全并行的寫回.最后的寫回請求的次數為group次，寫回時間穩定在2～3 s 且不隨核數的增加而增加，十分的穩定，不隨超算集群的用戶數量、網絡環境的波動而波動.

4 計算、讀取通信、寫回3 層重疊

4.1 讀取通信與計算重疊

由于在具體的數據同化過程中，每個計算進程只要獲取到相應的數據，同化過程便可以進行，不需要其他的通信操作，各個計算進程完全并行.針對這種同化過程，本文實現了對讀取通信和計算的重疊.為了實現讀取通信與計算的重疊，Reader 將需要讀取的數據進行緯度ny方向的切片，生成overlap片數據塊.每次Reader 只讀取1 個小塊的數據，在讀取到相應的數據后，再進行Reader 到Master、Master 到節點內計算進程的通信，各個計算進程在獲取到數據后便進行同化計算.在其他進程進行同化計算的同時，Reader 再讀取下一小片的數據，并和Master 進行通信，這樣便實現讀取通信和計算的重疊.

由于同化的數據都是帶邊的，其中經度nx方向的prep_rx＝105，緯度ny方向的prep_ry＝19，假設每個計算進程所需數據塊的nx方向子區域的長度為sub_nx_len，ny方向的區域分組的長度為ny_len，那么每個Reader 需要讀 取的數據為（（nx+2×prep_rx）×（ny_len+2×prep_ry））×sizeof（type），每個計算進程所需的通信數據大小為（sub_nx_len+2×prep_rx×（ny_len+2×prep_ry））×sizeof（type）.因 此，ny方 向上數據 切得越細，需要的帶邊數據也就越多，所產生的讀取和通信冗余量就越大.

為了解決讀取和通信數據冗余的問題，本文采用在讀取時先提前開辟（（nx+2×prep_rx）×（ny_len+2×prep _ry））×sizeof（type）的內存空間，其大小為不分片時的整體數據塊的內存大小加上帶邊數據的內存大小，其中ny_len為區域分組中數據塊在ny方向上的總長度.然后設置了psi_read_index作為Reader 已經讀取的數據偏移，在Reader 進行讀取時根據psi_read_index計算下一次需要讀取的數據量，即邏輯上需要的數據偏移減去psi_read_index代表的已經讀入內存的數據偏移.在實際執行時，分片后數據讀取所讀入的數據量和以不分片方式讀入的數據量相等，完全解決了帶邊數據導致的讀取數據冗余問題，實現了高效的數據復用.

Fig.13 The overlap between reading,communication and computation圖 13 讀取、通信與計算重疊

解決完讀取數據冗余的問題，下面介紹通信數據冗余問題的解決.本文算法在通信之前，和讀取相似，先提前開辟（sub_nx_len+2×prep_rx×（ny_ len+2×prep_ry））×sizeof（type）的內存空間，并設置psi_index為已經通過通信獲取到的同化數據的數據偏移.在通信時利用psi_index計算需要通信的數據量，即將需要通信的數據偏移減去已經通過通信的獲取到的數據偏移.通過這種方式，最終的通信數據量和未分片的數據量相等，完全解決了帶邊數據所導致的通信數據冗余問題.

具體的讀取、通信與計算重疊算法和避免讀取、通信冗余策略如圖13、圖14 所示.

4.2 計算、讀取通信、寫回重疊

2.1節中提到計算拓撲圖可以反映出本次數據同化計算任務的分布特性.同化數據顯示的是全球的海洋氣候信息，所以包含大量的陸地以及大量沒有觀測點的海洋區域，而且這些區域多是成塊且連續的.因此根據海洋同化數據計算量的分布特性，設計了讀取和通信避免策略.

在生成計算拓撲圖后，對整個海洋同化數據進行分類，分為需要同化區域與不需要同化區域.對于需要同化區域則根據上文的負載均衡算法進行基于負載均衡的任務劃分以及同化數據讀取、分發和計算.對于不需要同化的區域，則可以設置一些Reader進行讀取和寫回，這部分區域雖然不用計算，但是為了生成全球海洋氣候的分析場，這些數據也應該寫回到最后的分析場中.

Fig.14 The strategy of avoiding reading and communication redundancy圖 14 避免讀取、通信冗余策略

在本文的實驗中，不需要同化的區域占總海洋同化數據的60%左右.對于負責這部分讀取和寫回的Reader 與其他進程相互獨立，可以獨立執行.因此在生成計算量分布圖后，這部分進程可以在其他進程進行計算和讀取通信的同時，并行地對不需要同化的數據進行寫回.這樣的數據讀取和寫回，可以完全被其他進程的同化計算所掩蓋掉，可以減少整體程序大量的數據讀取通信和寫回的時間.經過上述的優化，計算進程由于只需要讀取通信需要同化的40%的區域，因此整體的讀取和通信量大大減少，讀取和通信時間也大大減少.與此同時，雖然有一些海洋區域不需要進行同化，但是如果將其分發給計算節點，由于計算節點需要對海洋格點逐個進行遍歷因此仍會存在一些判定處理的時間，在將這部分不需要同化計算的區域篩選出來后，只將需要同化計算的部分劃分給計算節點，計算節點整體的同化計算效率也大大提升，同化計算時間也大大減少.

5 實驗分析

本文采用的實驗平臺分別是天河2 號超算集群和曙光7000 超算集群，這2 個平臺所使用的數據集參數分別如表1 和表2 所示.

5.1 天河上各部分性能對比測試

1）讀取時間測試

本文分別對I/O 代理程序的讀取通信策略和基于Lustre 文件存儲架構的、基于3 層通信的讀取通信策略進行了對比測試.如圖15 所示，I/O 代理程序的讀取較為耗時，本文程序的讀取時間不僅遠小于I/O 代理程序的讀取時間，且隨核數的增加而基本保持在2 s 左右.

Table 1 Parameters of Data Set Used on Tianhe表 1 天河上使用的數據集參數

Table 2 Parameters of Data Set Used on Sugon 7000表 2 曙光7000 上使用的數據集參數

Fig.15 Comparison of reading time between the two algorithms圖 15 2 種算法讀取時間對比

2）通信時間測試

從圖16 可以看出I/O 代理程序的通信時間隨核數的增加而線性增加，而本文程序的通信時間隨著核數的增加基本保持不變.

Fig.16 Comparison of communication time between the two algorithms圖 16 2 種算法通信時間對比

3）數據同化時間測試

本文負載均衡優化前后的數據同化時間對比如圖17 所示.基于計算拓撲圖的負載均衡算法充分利用了可用的核數，避免了少部分核被分配大量計算量而大部分核只負責少部分甚至沒有分配計量的情況，充分利用了超算集群的計算性能，使得同化時間大大縮短.

Fig.17 Comparison of assimilation time between the two algorithms圖 17 2 種算法同化時間對比

4）寫回時間測試

本文分別對2 種算法的寫回時間做了測試，如表3 所示.I/O 代理程序的寫回時間十分的長，而且在實際測試時表現得十分不穩定，本文程序的寫回時間不僅遠小于I/O 代理程序的寫回時間且隨著核數的增加基本保持穩定.

Table 3 Comparison of Write Back Time Between the Two Algorithms表 3 2 種算法的寫回時間對比

5）預同化時間測試

本文預同化過程涉及在千核下5.4 GB 文件的讀取和分發.從圖18 可以看出.在未優化前，預同化時間隨著核數的增加呈線性增長；優化后隨著核數的增加預同化時間基本保持穩定.

Fig.18 Comparison of pre-assimilation time before and after optimization圖 18 優化前后預同化時間對比

6）整體性能測試

本文將展示整體并行同化程序優化前后的運行時間對比.如圖19 所示，分別顯示了I/O 代理程序、本文非計算重疊版本和本文最終優化版本的程序整體運行時間對比.從圖中可以看出I/O 代理程序的計算、讀取和寫回時間都比較大.雖然計算時間隨著核數的增加而明顯減少，但是文件的I/O 讀取通信時間隨著核數的增加而逐漸增加.同時I/O 代理程序的寫回時間隨著核數增加也十分不穩定.在非計算重疊算法中，可以看出程序的寫回時間和預同化時間仍占據一定的比重，且預同化時間隨著核數的增加而增加.在本文最終的優化版本中除了計算、讀取通信和寫回重疊的運行時間以外，其他部分的運行時間占總時間的比例不大.而且對比前后2 種優化算法，可以明顯看出讀取通信以及寫回的時間，在計算時間足夠長時實現了完全的重疊.同時也可以看出由于篩選出了60%不需要同化計算的區域，讀取和通信的數據量大大減少，本文最終優化版本的讀取和通信時間大大減少.而且在同化計算中由于省去了60%同化區域判定處理的時間，整體海洋數據同化程序的同化計算時間也大大減少.總體來說經過重疊優化后本文最終優化版本相比于之前優化的版本性能得到了很大的提升.在4 000 核下本文程序的整體性能相比I/O 代理程序提升了18 倍.

Fig.19 The overall running time of thousand cores program on Tianhe圖 19 天河上千核程序整體運行時間

5.2 曙光7000 上整體性能對比測試

Fig.20 The overall running time of thousand DCU card program on Sugon圖 20 曙光上千塊DCU 卡程序整體運行時間

在曙光7000 超算集群上本文算法仍然表現優異.由于節點布局與天河不同，本文在節點分配上做了調整.曙光7000 是Parastor 文件系統，并且未提供用戶處理節點（類似于Lustre 文件系統的OST）的接口，但節點（OST）、節點池（OSS）眾多，更加偏向于隨機散列的讀取，因此本文在文件讀取和通信優化上取消了Reader 這一層，采用Master 直接和OST 進行交互，以最大限度地利用Parastor 節點OST 的獨立通道，提高并行讀取效率.取消Reader 這一層，Master直接和Split 組進行交互，也將提升Master 到計算進程的通信效率.在圖20 中分別列出了預同化、讀取、通信、同化計算和寫回時間，其他比重過小的時間這里未做展示.由于曙光7000 超算集群上的數學庫還不完善，本文使用的MAGMA 數學庫計算效率很低，雖然使用了DCU 卡進行加速，但是其計算速度仍慢于天河上的MKL 庫函數，但這并不影響本文算法程序與I/O 代理程序性能的對比.在使用相同庫函數的情況下本文算法程序性能明顯優于I/O 代理程序，在4 000 塊卡下整體性能提升8 倍左右.

6 結束語

本文主要介紹了基于LICOM 模式的海洋數據同化程序的優化算法.分別提出了基于計算拓撲圖的負載均衡策略，基于Lustre 文件存儲架構的讀取、寫回和通信策略，并在此基礎上實現了計算、讀取通信、寫回的3 層重疊.相比于I/O 代理程序，本文程序的各項性能都得到了顯著的提升.在天河上，4 000 核下本文程序整體性能提升18 倍.在曙光7000 上，4 000塊卡下，整體性能提升了8 倍左右.

但是本程序仍有可優化的空間，由于隨著核數的增加計算時間不斷地縮短，而I/O 時間隨著核數的增加不會減少，因此當核數足夠大時計算的時間會小于I/O 時間，計算時間無法再掩蓋I/O 時間，這樣會限制程序的可拓展性上限.由于曙光7000 的存儲架構以及節點分布特性與天河不同，再加上曙光7000文件系統用戶接口的限制，雖然本文針對其做了相應的優化，盡可能地將Parastor 文件系統的節點（OST）利用起來，但是其仍達不到像天河Lustre 文件系統一樣極限的I/O 帶寬利用率.因此在未來的工作中準備進一步優化I/O 時間以提升程序的可拓展性上限，并對曙光7000 的文件存儲架構做進一步的研究，對曙光7000 上的程序做進一步的優化.

作者貢獻聲明：蔡迪負責論文相關的實驗設計、代碼實現、實驗測試以及論文撰寫等工作；洪學海負責前期實驗方案的設計討論與論文修改；肖俊敏負責論文結構的討論和修改；譚光明負責實驗補充和設計討論，并指導論文修改.