THREP：面向地球系統(tǒng)模式的插值算法研究平臺(tái)

吳 竑，王小鴿，楊廣文，朱昌磊，宋順強(qiáng)

清華大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)系，北京 100084

THREP：面向地球系統(tǒng)模式的插值算法研究平臺(tái)

吳 竑，王小鴿，楊廣文，朱昌磊，宋順強(qiáng)

清華大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)系，北京 100084

1 引言

1.1 背景

插值算法[1]是一類應(yīng)用十分廣泛的算法，它們?cè)趫D形圖像處理、數(shù)值計(jì)算等許多科學(xué)研究領(lǐng)域都是不可或缺的重要數(shù)學(xué)手段。插值算法在耦合地球系統(tǒng)模式中也有著極其重要的應(yīng)用。

耦合地球系統(tǒng)模式[2]是一套對(duì)地球系統(tǒng)進(jìn)行定量研究的軟件工具，地球科學(xué)相關(guān)領(lǐng)域的科學(xué)家根據(jù)專業(yè)知識(shí)對(duì)地球各個(gè)子系統(tǒng)進(jìn)行建模，編寫(xiě)各自領(lǐng)域的模擬程序，再由這些模擬程序相互耦合，成為一個(gè)完整的地球系統(tǒng)模式。為了使不同研究人員的代碼能夠耦合起來(lái)形成一個(gè)整體，往往需要他們共同遵循某個(gè)耦合軟件框架[3]。目前較為流行的地球系統(tǒng)模式耦合框架有歐盟PRISM所采用的OASIS[4]和美國(guó)CCSM3所采用的COUPLER[5]系列等。然而，由于地學(xué)領(lǐng)域的科學(xué)家們對(duì)地球各子系統(tǒng)進(jìn)行建模時(shí)所使用的數(shù)值方法不同，導(dǎo)致他們?cè)诰幊虝r(shí)使用不同的地球網(wǎng)格劃分方式（如經(jīng)緯網(wǎng)格、高斯網(wǎng)格[6]），從而對(duì)地球系統(tǒng)模式中不同分量模式間的數(shù)據(jù)交換提出了挑戰(zhàn)。因此，在不同網(wǎng)格間進(jìn)行數(shù)據(jù)交換時(shí)需要一種數(shù)據(jù)映射方式——插值。此時(shí)，可以將每個(gè)網(wǎng)格單元看做一個(gè)插值點(diǎn)，并用網(wǎng)格單元中心點(diǎn)坐標(biāo)表示插值點(diǎn)坐標(biāo)。

1.2 相關(guān)工作

SCRIP[7]是一個(gè)面向地球系統(tǒng)模式的插值算法庫(kù)，可以離線地生成地球系統(tǒng)模式研究人員所需地球網(wǎng)格間的插值權(quán)重系數(shù)文件，它被使用于CPL6等耦合器中。但是SCRIP中的插值算法不很完善，其守恒插值算法只支持一階，雙三次插值算法也沒(méi)有實(shí)現(xiàn)求導(dǎo)數(shù)的方案。此外，由于SCRIP中的數(shù)據(jù)結(jié)構(gòu)過(guò)于簡(jiǎn)單，它能夠支持的網(wǎng)格類型非常有限，這些因素都制約了其在地球系統(tǒng)模式中的作用。

ΤHU-Remap（https：//github.com/xunzhang/Τ-Remap）是針對(duì)SCRIP的缺陷，開(kāi)發(fā)的一套全新插值算法庫(kù)。它一方面重構(gòu)了SCRIP中一些算法，如反距離權(quán)重算法、雙線性算法，大大提高了它們的計(jì)算性能；另一方面改進(jìn)了SCRIP中的一些不足，如實(shí)現(xiàn)了雙三次差值中的求導(dǎo)、提供守恒插值中閾值的自動(dòng)調(diào)整等。但是，ΤHU-Remap對(duì)非規(guī)則網(wǎng)格的支持依然不夠，這不利于分量模式與耦合器的長(zhǎng)期發(fā)展。

ESMF（Earth System Modeling Framework）[8]項(xiàng)目由NASA地球科學(xué)辦公室于2000年9月提出，是一套面向地球系統(tǒng)模式的并行軟件開(kāi)發(fā)框架。它不含任何分量模式代碼，而是提供標(biāo)準(zhǔn)框架接口，用戶可以基于此接口實(shí)現(xiàn)天氣、氣候、數(shù)據(jù)同化等相關(guān)領(lǐng)域的應(yīng)用模型并將它們耦合起來(lái)成為一個(gè)完整的地球系統(tǒng)模式。ESMF模塊化設(shè)計(jì)使其具有極高的可重用性、互操作性和可移植性，這將有利于模式專家專注于模式本身的發(fā)展并減少其耦合開(kāi)銷(xiāo)。ESMF主要提供Superstructure和Infrastructure兩部分功能供開(kāi)發(fā)者使用：Superstructure中封裝了地球系統(tǒng)模式中的耦合驅(qū)動(dòng)和接口，也模塊化了分量模式中物理過(guò)程、動(dòng)力框架及一些子過(guò)程；Infrastructure則提供更加底層的抽象，它又可細(xì)分為數(shù)據(jù)處理和公用模塊兩部分，數(shù)據(jù)處理部分主要封裝所有與地球系統(tǒng)模式中數(shù)據(jù)處理相關(guān)的數(shù)據(jù)結(jié)構(gòu)及功能，如Array、Grid、Field等，公用模塊部分則提供了如數(shù)值計(jì)算、時(shí)間管理、進(jìn)程管理、I/O、容錯(cuò)日志處理的基本功能。

ESMF-Regridding[9]是ESMF提供的一個(gè)基于ESMF框架開(kāi)發(fā)的插值應(yīng)用，它能并行生成離線插值權(quán)重系數(shù)文件，特別是它能支持非規(guī)則網(wǎng)格間的插值。但是，ESMF-Regridding僅以一個(gè)小工具的形式出現(xiàn)在ESMF中，沒(méi)有抽象插值在應(yīng)用層面的功能供開(kāi)發(fā)人員使用，因此不夠靈活和通用。插值算法開(kāi)發(fā)者往往需要閱讀大量ESMF源碼，作出相應(yīng)修改并實(shí)現(xiàn)各自新的插值算法，難度極大，違背了ESMF的設(shè)計(jì)目標(biāo)。此外，ESMF-Regridding沒(méi)有提供驗(yàn)證插值算法正確性與性能的功能，不利于深入研究與對(duì)比各插值算法的特性。

雖然目前已有SCRIP、ΤHU-Remap、ESMF-Regridding等面向地球系統(tǒng)模式的插值軟件，但隨著地球系統(tǒng)模式研究的深入，對(duì)插值算法提出了新要求。為插值算法的研究、測(cè)試和評(píng)價(jià)建立一個(gè)統(tǒng)一的平臺(tái)，具有重要的實(shí)踐意義。基于ESMF框架，開(kāi)發(fā)了一個(gè)通用、靈活、并行的插值研究平臺(tái)ΤHREP（Τsinghua Regridding Platform），它能同時(shí)滿足插值算法研究人員、耦合器開(kāi)發(fā)者、模式專家三類用戶需求。ΤHREP集成了ΤHU-Remap與ESMF-Regridding中的插值算法，封裝了易用的插值接口，并提供了驗(yàn)證、對(duì)比等功能。

2 平臺(tái)設(shè)計(jì)

ΤHREP綜合考慮了算法研究人員、耦合器開(kāi)發(fā)者及模式專家三類人對(duì)插值的需求，是一個(gè)通用、靈活、并行的面向地球系統(tǒng)模式的插值平臺(tái)。它抽象了插值過(guò)程，封裝了插值流程，便于已有插值算法庫(kù)的接入，利于長(zhǎng)期發(fā)展和標(biāo)準(zhǔn)化。同時(shí)，為滿足不同用戶，ΤHREP在保證易用性的同時(shí)，從底到上提供了多層接口。ΤHREP還提供了插值檢驗(yàn)功能，便于深入研究各算法的精度、性能等相關(guān)問(wèn)題。

2.1 平臺(tái)架構(gòu)

ΤHREP插值平臺(tái)的總體結(jié)構(gòu)圖見(jiàn)圖1，主要包括輸入文件、配置模塊、權(quán)重生成模塊、插值模塊和檢驗(yàn)?zāi)K幾個(gè)部分。輸入文件含網(wǎng)格數(shù)據(jù)文件和物理量數(shù)據(jù)文件，配置模塊負(fù)責(zé)配置一次插值過(guò)程需要的信息，權(quán)重生成模塊包含了初始化和搜索兩個(gè)子模塊，ΤHREP通過(guò)初始化模塊建立相應(yīng)數(shù)據(jù)結(jié)構(gòu)的對(duì)象，在此基礎(chǔ)上搜索找出各算法所需格點(diǎn)信息，生成插值權(quán)重，插值模塊通過(guò)矩陣向量乘計(jì)算出實(shí)際物理量插值后的結(jié)果從而完成整個(gè)插值過(guò)程，檢驗(yàn)?zāi)K將對(duì)插值的精度和性能作進(jìn)一步分析并畫(huà)出圖形。

圖1 ΤHREP平臺(tái)架構(gòu)

2.2 通用性

目前地球系統(tǒng)模式中常用的插值算法有反距離權(quán)重插值[7]、雙線性插值、雙三次插值[10]、Patch插值[11-12]、三次樣條插值[13]、守恒插值[7，14]等，ΤHREP里將所有插值算法抽象成了如下兩個(gè)步驟：

（1）計(jì)算插值權(quán)重系數(shù)矩陣。這個(gè)階段又可細(xì)分為兩個(gè)子過(guò)程：首先根據(jù)源網(wǎng)格和目標(biāo)網(wǎng)格間的拓?fù)潢P(guān)系，對(duì)每個(gè)目標(biāo)網(wǎng)格單元（cell）找到其對(duì)應(yīng)位置在源網(wǎng)格的“若干”網(wǎng)格單元（這里的“若干”由不同插值算法決定，如雙線性插值要找出包住目標(biāo)網(wǎng)格單元的4個(gè)源網(wǎng)格單元，即目標(biāo)網(wǎng)格在找出的4個(gè)源網(wǎng)格單元連線構(gòu)成的四邊形內(nèi)部）。然后根據(jù)插值算法計(jì)算出每個(gè)找到的源網(wǎng)格單元相對(duì)于該目標(biāo)網(wǎng)格單元的權(quán)重。因?yàn)槊總€(gè)目標(biāo)網(wǎng)格單元都要計(jì)算其對(duì)應(yīng)的“若干”個(gè)權(quán)重值，因此這一步完成后，會(huì)得到一個(gè)權(quán)重矩陣W，W的每行對(duì)應(yīng)每個(gè)目標(biāo)網(wǎng)格單元id，每列對(duì)應(yīng)每個(gè)源網(wǎng)格單元id。這里W是一個(gè)稀疏矩陣，如雙線性插值的權(quán)重矩陣每行只有4個(gè)非零元素。

ΤHREP中使用NetCDF[15]標(biāo)準(zhǔn)的文件格式作為整個(gè)平臺(tái)數(shù)據(jù)流中輸入、輸出文件的格式。NetCDF文件是一種格式化的二進(jìn)制文件，它既有二進(jìn)制文件的數(shù)據(jù)連續(xù)性，又能使用戶像文本文件那樣方便地讀取相應(yīng)屬性的數(shù)據(jù)。采用NetCDF文件的另一個(gè)重要原因是它已經(jīng)是研究地球系統(tǒng)模式人員的約定標(biāo)準(zhǔn)，分量模式中使用的網(wǎng)格文件均采用這種格式，因此用NetCDF作為插值結(jié)果的輸出文件能更便于用戶直接使用。

考慮通用性的另一方面是功能可擴(kuò)展，而數(shù)據(jù)結(jié)構(gòu)直接決定了能否支持眾多功能。因此，通用的數(shù)據(jù)結(jié)構(gòu)也是系統(tǒng)通用性的關(guān)鍵。整個(gè)插值過(guò)程中主要的數(shù)據(jù)包括網(wǎng)格信息和物理量信息，插值所用的數(shù)據(jù)結(jié)構(gòu)主要就是為了便于插值計(jì)算而在這兩部分?jǐn)?shù)據(jù)信息基礎(chǔ)上建立的。在ΤHREP中盡量保證了數(shù)據(jù)的原始邏輯，并刻畫(huà)了網(wǎng)格單元間的相鄰邏輯關(guān)系，這樣可以保證功能可擴(kuò)展，便于平臺(tái)的進(jìn)一步研發(fā)，也為以后我們或者其他開(kāi)發(fā)者封裝更高級(jí)的數(shù)據(jù)結(jié)構(gòu)提供了可能。

2.3 靈活性

在眾多地球系統(tǒng)模式的研發(fā)和使用人員中主要有三類人員涉及插值。第一類是模式專家，他們負(fù)責(zé)開(kāi)發(fā)分量模式，專注于各分量模式的功能、精度與性能問(wèn)題。因各分量模式在自然界中并不是孤立存在，因此一個(gè)模式依賴另一個(gè)模式中插值得來(lái)的數(shù)據(jù)做計(jì)算；第二類是耦合器開(kāi)發(fā)人員，耦合器負(fù)責(zé)分量模式之間的數(shù)據(jù)管理，而數(shù)據(jù)管理主要包括通信、插值、集合等子模塊，因此耦合器開(kāi)發(fā)人員也關(guān)注插值；第三類是插值算法研究人員，他們專注于插值算法本身的研究，如SCRIP就是他們?yōu)榈厍蛳到y(tǒng)模式中的離線插值開(kāi)發(fā)的一個(gè)算法庫(kù)。為了能同時(shí)滿足三類研究人員的需求，ΤHREP在對(duì)插值流程進(jìn)行封裝的基礎(chǔ)上提供了并行生成權(quán)重矩陣、計(jì)算插值、驗(yàn)證插值正確性、插值性能測(cè)試等功能。

生成權(quán)重矩陣是離線生成一組網(wǎng)格對(duì)應(yīng)的插值系數(shù)文件，對(duì)于一組網(wǎng)格（一個(gè)源網(wǎng)格和一個(gè)目標(biāo)網(wǎng)格）只需計(jì)算一次。在線耦合過(guò)程將調(diào)用該系數(shù)文件，再完成插值的第二步矩陣向量乘操作，而這個(gè)系數(shù)文件可能被多次使用。同時(shí)，ΤHREP支持直接計(jì)算插值，即不產(chǎn)生權(quán)重系數(shù)文件，直接進(jìn)行插值計(jì)算，得到插值結(jié)果，ΤHREP中計(jì)算插值權(quán)重的過(guò)程是并行的。平臺(tái)還提供了插值性能測(cè)試功能，為研究高性能插值提供了可能。最后，驗(yàn)證模塊主要基于插值結(jié)果文件，用兩套驗(yàn)證策略對(duì)計(jì)算出的插值結(jié)果進(jìn)行精度分析，并計(jì)算一些統(tǒng)計(jì)量信息、同時(shí)畫(huà)出相應(yīng)的誤差分布圖。

隨著分量模式自身和整個(gè)地球系統(tǒng)模式的發(fā)展，在線插值越來(lái)越被重視。之所以要把離線的權(quán)重計(jì)算操作放到在線去做，是因?yàn)楦鞣至磕Ｊ降木W(wǎng)格在模式迭代計(jì)算過(guò)程中可能不斷發(fā)生變化，這樣就要頻繁地計(jì)算插值權(quán)重系數(shù)。這時(shí)對(duì)耦合器來(lái)說(shuō)，插值的兩步之間不再需要產(chǎn)生權(quán)重系數(shù)文件，而是直接算出插值結(jié)果，ΤHREP目前只支持串行這么做。當(dāng)然，在線插值對(duì)插值性能提出了很高的要求，如何并行支持在線插值有待進(jìn)一步研究。

ΤHREP還具有框架的特性，用戶可以方便地在其接口上實(shí)現(xiàn)新插值算法。同時(shí)，ΤHREP中數(shù)據(jù)結(jié)構(gòu)的通用性使開(kāi)發(fā)插值算法庫(kù)的人能方便地作轉(zhuǎn)換，將已有插值庫(kù)快速接入到ΤHREP平臺(tái)中。ΤHREP也提供了插值算法層的接口，用戶能方便快速地在平臺(tái)中開(kāi)發(fā)新的插值算法。ΤHREP中將插值分為結(jié)構(gòu)網(wǎng)格間的插值與非結(jié)構(gòu)網(wǎng)格間的插值，由于結(jié)構(gòu)網(wǎng)格的索引能更清晰地反映網(wǎng)格的拓?fù)潢P(guān)系，為使用簡(jiǎn)單，將結(jié)構(gòu)網(wǎng)格的接口與非結(jié)構(gòu)網(wǎng)格的接口分離開(kāi)來(lái)，而底層實(shí)現(xiàn)則采用一致的數(shù)據(jù)結(jié)構(gòu)。

2.4 并行性

并行生成插值權(quán)重系數(shù)是ΤHREP的另一特點(diǎn)，這為研究插值算法性能的用戶提供了有利工具，也為以后支持在線插值提供了可能?？紤]到代碼復(fù)用，不想讓插值算法研究人員為了實(shí)現(xiàn)并行做重構(gòu)，因此希望并行盡量地普適。并行生成插值權(quán)重本質(zhì)上是對(duì)源網(wǎng)格與目標(biāo)網(wǎng)格進(jìn)行剖分，這對(duì)于分布式環(huán)境中的每個(gè)獨(dú)立計(jì)算單元P來(lái)看，計(jì)算方式和串行情況下幾乎一致，而區(qū)別在于此時(shí)的網(wǎng)格不再是整個(gè)地球展開(kāi)的平面，而是剖分后的一部分區(qū)域。如果把生成插值權(quán)重系數(shù)的過(guò)程抽象成一個(gè)算子Gen，那么唯一造成需要Gen重構(gòu)的可能性是Gen與網(wǎng)格的拓?fù)湫螤钣嘘P(guān)。為此，將網(wǎng)格做了規(guī)則的矩形剖分，從而保持原有的網(wǎng)格形狀，由用戶來(lái)負(fù)責(zé)實(shí)現(xiàn)通信邏輯。ΤHREP提供一維通信的接口，用戶可以使用它簡(jiǎn)化通信細(xì)節(jié)的實(shí)現(xiàn)。

2.5 插值檢驗(yàn)

作為研究平臺(tái)，ΤHREP具有插值檢驗(yàn)功能。為了驗(yàn)證精度，需要對(duì)插值結(jié)果的正確性作定量分析?，F(xiàn)設(shè)νalcal為ΤHREP中插值得出的結(jié)果，νalacc為物理量數(shù)據(jù)的精確值，則目標(biāo)插值格點(diǎn)插值結(jié)果的相對(duì)誤差er可用公式er= |νalacc-νalcal|/νalacc算出。

問(wèn)題的難點(diǎn)在于如何刻畫(huà)公式中的νalacc，因?yàn)閷?shí)際上并不知道目標(biāo)插值網(wǎng)格上的真實(shí)物理量分布數(shù)據(jù)，為此ΤHREP中采用了兩種方式來(lái)定義νalacc。一是用解析函數(shù)值，即給定解析函數(shù)表達(dá)式，要求某坐標(biāo)處的νalacc只需直接取函數(shù)值便可。這時(shí)為了盡可能真實(shí)地模擬不同物理量在不同地形的分布特性，可以取球諧函數(shù)[16]作為驗(yàn)證的解析函數(shù)。然而，真實(shí)的物理量分布式不能簡(jiǎn)單用解析函數(shù)刻畫(huà)的，用解析函數(shù)測(cè)試精度較高的插值算法在真實(shí)物理量分布中不一定仍然精確。因此ΤHREP還使用另一種策略來(lái)定義νalacc：先將真實(shí)物理量的月平均數(shù)據(jù)用守恒插值插到一個(gè)“足夠密”的網(wǎng)格上（“足夠密”的含義是對(duì)于任何插值網(wǎng)格上的每一個(gè)點(diǎn)，都存在一個(gè)該“足夠密”網(wǎng)格上的格點(diǎn)，使得它們之前的距離“足夠近”），忽略這次插值的誤差。這樣，就可以用此“足夠密”網(wǎng)格格點(diǎn)上的值代表目標(biāo)插值網(wǎng)格上的物理量的值了。ΤHREP中定義的“足夠密”網(wǎng)格分辨率為1 800×900，因此上述“足夠近”大約在3.5 km范圍內(nèi)。這兩種確定νalacc的方式分別對(duì)應(yīng)著兩套檢驗(yàn)策略，用戶在進(jìn)行完插值后，可以選擇任意一種進(jìn)行檢驗(yàn)，平臺(tái)將通過(guò)NCL[17]腳本畫(huà)出誤差分布圖，并同時(shí)算出平均相對(duì)誤差（ΑVG_ERR）、最大相對(duì)誤差（MΑX_ERR）、均方根誤差（RMSE）等統(tǒng)計(jì)信息。

2.6 數(shù)據(jù)流

ΤHREP的數(shù)據(jù)流圖如圖2所示。圖2（a）是頂層圖，首先ΤHREP通過(guò)讀入配置參數(shù)和網(wǎng)格信息產(chǎn)生插值算子對(duì)象和網(wǎng)格對(duì)象，同時(shí)讀入源網(wǎng)格上物理量信息，平臺(tái)利用這些對(duì)象與信息進(jìn)行插值，從而求出插值結(jié)果。驗(yàn)證模塊讀入插值結(jié)果，畫(huà)出物理量分布圖與插值誤差分布圖。

圖2（a）ΤHREP數(shù)據(jù)流頂層圖

圖2（b）ΤHREP數(shù)據(jù)流0層圖

在整個(gè)插值計(jì)算模塊內(nèi)部，數(shù)據(jù)流向如0層圖2（b）所示，數(shù)據(jù)將依次流過(guò)兩個(gè)子模塊：權(quán)重生成模塊與計(jì)算插值模塊。權(quán)重生成模塊負(fù)責(zé)完成插值權(quán)重的計(jì)算，若配置要求離線插值，則此步將產(chǎn)生一個(gè)權(quán)重系數(shù)文件供耦合開(kāi)發(fā)者戶和模式專家使用，且此過(guò)程是并行的，若配置要求的是在線插值，ΤHREP則將在完成權(quán)重計(jì)算后直接計(jì)算出插值，得到物理量在目標(biāo)網(wǎng)格單元上的值。需要指出的是，計(jì)算插值這個(gè)子過(guò)程本質(zhì)上是一個(gè)稀疏矩陣乘向量的過(guò)程，雖然它本身易并行，但這時(shí)的并行剖分并不同于最初讀入網(wǎng)格的剖分，而是需要進(jìn)行一次全局通信來(lái)產(chǎn)生新剖分。因此，若想對(duì)整個(gè)插值計(jì)算過(guò)程做并行，則務(wù)必需要進(jìn)行一次全局通信，隨著剖分?jǐn)?shù)變大，并行的效率也將大大下降?？紤]到效率，ΤHREP中并沒(méi)有提供并行地支持在線插值而只是并行地生成插值權(quán)重系數(shù)，要解決在線插值的并行效率，必須使得源網(wǎng)格的剖分充分靈活，從而減少甚至消除內(nèi)部的全局通信。另外，對(duì)于插值結(jié)果，如圖2（b）所示，驗(yàn)證模塊可以通過(guò)解析函數(shù)和真實(shí)物理量?jī)煞N定義νalacc的策略分別進(jìn)行驗(yàn)證。

3 平臺(tái)實(shí)現(xiàn)

地球系統(tǒng)模式中的插值不同于圖形學(xué)中插值，由于它涉及的坐標(biāo)點(diǎn)都在地球表面，因此可以認(rèn)為是一種近似的平面插值。之所以說(shuō)近似是因?yàn)榈厍虮砻娌⒎瞧矫妫以跇O點(diǎn)和周期邊界處的處理也不同于普遍意義下的二維插值。本文基于上述設(shè)計(jì)，用C++實(shí)現(xiàn)了ΤHREP，并封裝了C++的插值接口供用戶使用。ΤHREP使用SCRIP標(biāo)準(zhǔn)網(wǎng)格文件格式，并以二維經(jīng)緯坐標(biāo)的形式存儲(chǔ)網(wǎng)格單元。為同時(shí)滿足上文中三類用戶的需求，ΤHREP采用了ESMF中Infrastructure的Mesh作為其內(nèi)部網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)，它足夠靈活和通用，符合ΤHREP的設(shè)計(jì)目標(biāo)。Mesh的數(shù)據(jù)結(jié)構(gòu)既能用于封裝結(jié)構(gòu)網(wǎng)格插值接口，又能用于實(shí)現(xiàn)非結(jié)構(gòu)網(wǎng)格間的插值。事實(shí)上，一個(gè)Mesh對(duì)象由nodes和elements構(gòu)成，每個(gè)node是一個(gè)帶坐標(biāo)信息的點(diǎn)，每個(gè)element代表由若干nodes構(gòu)成的一個(gè)單元區(qū)域。Element還給出了Mesh的形狀信息與nodes間的相鄰?fù)負(fù)潢P(guān)系。對(duì)于結(jié)構(gòu)網(wǎng)格，node就可以滿足插值需求，而對(duì)非結(jié)構(gòu)網(wǎng)格，則需加上element的信息。對(duì)于非規(guī)則網(wǎng)格插值的支持也意味著ΤHREP可用于帶mask的網(wǎng)格，如海洋網(wǎng)格或海冰網(wǎng)格。考慮到球面插值的特殊性，ΤHREP在網(wǎng)格坐標(biāo)對(duì)象中加入了額外的ghost點(diǎn)來(lái)處理周期邊界問(wèn)題，而對(duì)于極點(diǎn)問(wèn)題ΤHREP提供三個(gè)選項(xiàng)供用戶進(jìn)行選擇。

ΤHREP的實(shí)現(xiàn)層級(jí)圖如圖3所示，底層采用NetCDF格式作為輸入文件格式，之上使用ESMF中Infrastructure的Mesh數(shù)據(jù)結(jié)構(gòu)表示網(wǎng)格信息以保證通用，在Mesh層以上，ΤHREP平臺(tái)封裝了兩套接口分別用于結(jié)構(gòu)網(wǎng)格的插值（Grid接口）與非結(jié)構(gòu)網(wǎng)格的插值（Mesh接口）?；贕rid，ΤHREP集成了ΤHU-Remap，其中包括反距離權(quán)重插值算法、三次樣條插值算法和雙三次插值算法?；贛esh，ΤHREP集成了ESMF-Regridding，其中包括雙線性插值算法、一階守恒及二階守恒插值算法。插值算法開(kāi)發(fā)者也可以針對(duì)其算法特性，選取Mesh或Grid接口進(jìn)行新算法的開(kāi)發(fā)。對(duì)于那些有獨(dú)特?cái)?shù)據(jù)結(jié)構(gòu)要求的算法或算法庫(kù)，用戶可以直接用Mesh數(shù)據(jù)結(jié)構(gòu)作轉(zhuǎn)換，而后集成到ΤHREP平臺(tái)中。此外，ΤHREP中的驗(yàn)證模塊底層實(shí)現(xiàn)則是基于NCL的。

圖3 ΤHREP實(shí)現(xiàn)層級(jí)圖

3.1 網(wǎng)格文件

ΤHREP選取SCRIP格式（一種結(jié)構(gòu)化NetCDF格式）作為網(wǎng)格文件，其文件頭如下：

Dimensions屬性中的grid_size表示該網(wǎng)格的網(wǎng)格單元總數(shù)，grid_rank表示維度（grid_rank為2表示2D邏輯矩形網(wǎng)格，為1表示非結(jié)構(gòu)網(wǎng)格）。Variables屬性中包含grid_dim，它代表每一維的網(wǎng)格單元數(shù)，如對(duì)于Τ42網(wǎng)格就是[128×64]。grid_center_lat和grid_center_lon分別表示網(wǎng)格單元中心點(diǎn)的經(jīng)緯度坐標(biāo)，grid_corner_lat和grid_corner_lon分表表示網(wǎng)格單元的頂點(diǎn)坐標(biāo)。需要注意的是，一個(gè)網(wǎng)格單元的頂點(diǎn)數(shù)由grid_corners表示，除了守恒插值其他的插值算法均用中心點(diǎn)的坐標(biāo)代表一個(gè)網(wǎng)格單元。grid_imask是一個(gè)代表是否為海洋的整數(shù)，它決定了該網(wǎng)格單元是否需要參與插值計(jì)算，特別當(dāng)grid_imask為0時(shí)，該網(wǎng)格單元不參與插值計(jì)算。網(wǎng)格的坐標(biāo)單元可以是角度或弧度，讀入網(wǎng)格信息后，ΤHREP根據(jù)文件中坐標(biāo)單位標(biāo)識(shí)統(tǒng)一將坐標(biāo)信息轉(zhuǎn)換成角度進(jìn)行處理。

3.2 極區(qū)問(wèn)題與周期邊界

平面插值與球面插值有著本質(zhì)的區(qū)別，因此將球面插值近似成平面插值需要對(duì)球面距離、極點(diǎn)以及周期邊界作特殊處理。

首先，ΤHREP用兩點(diǎn)的直線距離近似兩點(diǎn)的球面大圓弧距離。一方面，大圓弧距離在計(jì)算量上遠(yuǎn)大于直線距離，另一方面，隨著分量模式發(fā)展，網(wǎng)格分辨率將越來(lái)越高，球面兩點(diǎn)的大圓弧距離與直線距離間的誤差將越來(lái)越小，插值結(jié)果也將越精確。

其次，對(duì)于極點(diǎn)處的插值，ΤHREP提供了三個(gè)選項(xiàng)：一是不處理，即遇到目標(biāo)網(wǎng)格單元中心點(diǎn)在源網(wǎng)格最低或最高緯度之外的情況時(shí)，運(yùn)行過(guò)程中將直接拋出異常。二是按正常情況處理而不拋出異常。這時(shí)需要插值算法自身來(lái)處理極點(diǎn)問(wèn)題，如現(xiàn)在平臺(tái)中的雙線性算法在處理目標(biāo)網(wǎng)格極點(diǎn)插值時(shí)，由于找不到包含它的源網(wǎng)格單元構(gòu)成的矩形，將會(huì)調(diào)用反距離權(quán)重的方式計(jì)算權(quán)重系數(shù)?？紤]到雙線性算法權(quán)重矩陣的大小，ΤHREP這時(shí)會(huì)找4個(gè)源網(wǎng)格單元并計(jì)算對(duì)應(yīng)權(quán)重。三是取所有包含極點(diǎn)的源網(wǎng)格單元值的平均值作為插值的結(jié)果。

最后，對(duì)于周期邊界，ΤHREP在生成網(wǎng)格數(shù)據(jù)對(duì)象的時(shí)候添加了ghost區(qū)域。所謂ghost區(qū)域就是經(jīng)度為0附近坐標(biāo)網(wǎng)格單元的集合，如經(jīng)度為359.5的網(wǎng)格點(diǎn)對(duì)應(yīng)的ghost坐標(biāo)點(diǎn)經(jīng)度為－0.5，而經(jīng)度為0.5的網(wǎng)格點(diǎn)對(duì)應(yīng)的ghost坐標(biāo)點(diǎn)經(jīng)度為360.5。事實(shí)上，這些ghost區(qū)域的點(diǎn)僅僅服務(wù)于插值所需的相鄰關(guān)系，而對(duì)于ghost格點(diǎn)單元上物理量數(shù)據(jù)，ΤHREP則將還原到原始網(wǎng)格對(duì)應(yīng)的網(wǎng)格單元上。

3.3 帶mask的插值

由于分量模式中使用的網(wǎng)格一般是對(duì)全球區(qū)域的剖分，如海洋模式中剖分后有些網(wǎng)格單元中將會(huì)包含陸地。這時(shí)，mask用來(lái)標(biāo)示該網(wǎng)格單元是否是海洋區(qū)域，因此對(duì)插值來(lái)說(shuō)，無(wú)論是源網(wǎng)格還是目標(biāo)網(wǎng)格，只要帶mask且 mask為1表示參與計(jì)算，而0表示不參與插值計(jì)算。

有些插值算法的實(shí)現(xiàn)依賴于網(wǎng)格拓?fù)湫螤?，如ΤHREP中的三次樣條插值算法，目前只支持矩形經(jīng)緯網(wǎng)格，因此一旦網(wǎng)格帶了mask，整個(gè)網(wǎng)格將不再是經(jīng)緯網(wǎng)格了。另外有些算法，ΤHREP中對(duì)源網(wǎng)格帶mask的情況作了特殊處理。如雙線性插值算法，其插值過(guò)程中需找出包含目標(biāo)網(wǎng)格單元的4個(gè)源網(wǎng)格單元，若此時(shí)找出的源網(wǎng)格單元被mask掉了，就不能使用，ΤHREP會(huì)調(diào)用反距離權(quán)重進(jìn)行插值，因?yàn)槠鋵?shí)現(xiàn)不受限于網(wǎng)格拓?fù)湫螤睢?/p>

3.4 并行

出于通用性考慮，ΤHREP中對(duì)源網(wǎng)格與目標(biāo)網(wǎng)格作的剖分均是不重疊的，這也就需要各插值算法進(jìn)一步實(shí)現(xiàn)通信細(xì)節(jié)。ΤHREP中的源網(wǎng)格剖分與目標(biāo)網(wǎng)格剖分方式是相互獨(dú)立的，目標(biāo)網(wǎng)格支持一維、二維剖分而源網(wǎng)格只支持一維剖分。對(duì)于距離權(quán)重算法與雙線性插值，由于計(jì)算每個(gè)目標(biāo)網(wǎng)格單元對(duì)應(yīng)的權(quán)重與源網(wǎng)格大小無(wú)關(guān)，而只要能“包含”目標(biāo)網(wǎng)格單元即可。因此只需對(duì)目標(biāo)網(wǎng)格進(jìn)行剖分，而可以不將源網(wǎng)格進(jìn)行剖分。這并不影響可擴(kuò)展性，估算分辨率為1°的經(jīng)緯網(wǎng)格（360×180）約只占0.5 MB的內(nèi)存。對(duì)于三次樣條等其他一些插值算法，需要同時(shí)對(duì)源網(wǎng)格與目標(biāo)網(wǎng)格進(jìn)行剖分。目前ΤHREP中的三次樣條插值支持對(duì)源網(wǎng)格進(jìn)行一維緯向剖分，對(duì)目標(biāo)網(wǎng)格進(jìn)行二維剖分。在插值開(kāi)始前，算法會(huì)將相鄰邏輯網(wǎng)格單元上的數(shù)據(jù)按層進(jìn)行通信，從而滿足上述“包含”的要求。

3.5 插值檢驗(yàn)?zāi)K

ΤHREP不僅能做插值計(jì)算，還能檢驗(yàn)插值結(jié)果的精度以及插值性能，這為研究插值提供了極大的便利。ΤHREP集成了NCL，使得能夠清晰直觀地畫(huà)出插值誤差分布圖。

有的插值算法在插值時(shí)不僅需要源網(wǎng)格上的函數(shù)值信息，還需要導(dǎo)數(shù)值信息。如二階守恒、雙三次插值等。本質(zhì)上，加入導(dǎo)數(shù)信息是為了更準(zhǔn)確刻畫(huà)真實(shí)物理量的變化趨勢(shì)。為測(cè)試這類插值算法，對(duì)于解析函數(shù)，只需直接算出導(dǎo)數(shù)表達(dá)式取對(duì)應(yīng)位置的導(dǎo)數(shù)值即可。而對(duì)于真實(shí)數(shù)據(jù)的方案，需要用數(shù)值的方法將物理量函數(shù)值近似出相應(yīng)的導(dǎo)數(shù)值。ΤHREP中采用了有限差分[18]和最小二乘擬合[19]的方法來(lái)近似導(dǎo)數(shù)。

總地來(lái)看，對(duì)于一次插值，ΤHREP能相應(yīng)地畫(huà)出插值分布圖、精確物理量的值分布圖與相對(duì)誤差分布圖。同時(shí)，ΤHREP也會(huì)計(jì)算出最大相對(duì)誤差（MΑX_ERR）、平均相對(duì)誤差（ΑVG_ERR）與均方根誤差（RMSE）等統(tǒng)計(jì)量信息。

3.6 插值接口與使用

ΤHREP還提供了方便的插值接口供用戶使用。為算法開(kāi)發(fā)人員提供了兩套接口，分別適用結(jié)構(gòu)網(wǎng)格與非結(jié)構(gòu)網(wǎng)格間的插值。對(duì)于結(jié)構(gòu)網(wǎng)格的插值，用戶只需繼承Common_remap基類，并在子類中實(shí)現(xiàn)兩個(gè)虛函數(shù)init_remap和cal_remap。同時(shí)，為了使添加的插值算法能并行執(zhí)行，ΤHREP封裝了幾種通信方式供用戶調(diào)用。ΤHREP會(huì)在整個(gè)插值過(guò)程中先后調(diào)用init_remap和cal_remap這兩個(gè)虛函數(shù)完成權(quán)重計(jì)算和插值計(jì)算。對(duì)于非結(jié)構(gòu)網(wǎng)格的插值，平臺(tái)提供了另一套接口幫助用戶滿足更復(fù)雜的插值需求。這套接口的通用性也為現(xiàn)有插值庫(kù)提供了接入的可能。事實(shí)上，ΤHREP中兩套算法接口底層是通過(guò)統(tǒng)一的數(shù)據(jù)結(jié)構(gòu)抽象出來(lái)的。ΤHREP中的算法注冊(cè)流程也是非常簡(jiǎn)單，便于新算法添加。最后，ΤHREP還為耦合器開(kāi)發(fā)者和模式專家提供了方便的執(zhí)行入口，只需簡(jiǎn)單的配置即可運(yùn)行使用。

3.6.1 結(jié)構(gòu)網(wǎng)格的插值接口

所謂的結(jié)構(gòu)網(wǎng)格是指其索引是M×N的，即每行N個(gè)網(wǎng)格單元，共M行。因此結(jié)構(gòu)網(wǎng)格的索引就能清晰地反映網(wǎng)格的拓?fù)潢P(guān)系，特別是網(wǎng)格單元的相鄰關(guān)系。針對(duì)這類網(wǎng)格的插值需求，提供了一個(gè)基類Common_remap如下：

用戶定義的類在初始化時(shí)設(shè)置set_size，Common_remap從而確定權(quán)重系數(shù)矩陣以及其他私有成員大小。定義的類必須繼承Common_remap類并實(shí)現(xiàn)init_remap和cal_remap兩個(gè)虛函數(shù)。

init_remap函數(shù)負(fù)責(zé)生成每個(gè)目標(biāo)網(wǎng)格單元的權(quán)重系數(shù)，其中需要調(diào)用write_wgts函數(shù)寫(xiě)入每個(gè)目標(biāo)網(wǎng)格單元對(duì)應(yīng)的set_size個(gè)權(quán)重，還需要調(diào)用write_matchup寫(xiě)入每個(gè)目標(biāo)網(wǎng)格單元對(duì)應(yīng)的set_size個(gè)源網(wǎng)格點(diǎn)信息。如雙線性插值，每個(gè)目標(biāo)網(wǎng)格單元產(chǎn)生4個(gè)源網(wǎng)格單元對(duì)應(yīng)的權(quán)重和4個(gè)源網(wǎng)格單元信息。Common_remap類還提供了兩個(gè)簡(jiǎn)化搜索的函數(shù)接口，get_nearest_nbrs能夠找出在threshold閾值范圍內(nèi)的若干個(gè)距離目標(biāo)網(wǎng)格單元最近的源網(wǎng)格單元。如距離權(quán)重插值、9個(gè)點(diǎn)的雙二次插值、16個(gè)點(diǎn)的雙三次插值等算法都能方便使用。而get_bnd_quadrangle能夠找出能夠包住目標(biāo)網(wǎng)格單元的由源網(wǎng)格單元構(gòu)成的四邊形，供雙線性插值、4個(gè)點(diǎn)的雙三次插值等一類算法使用。

cal_remap負(fù)責(zé)通過(guò)init_remap算出的權(quán)重矩陣與源網(wǎng)格單元上物理量值做插值。事實(shí)上，對(duì)于大部分算法，該過(guò)程是一個(gè)簡(jiǎn)單的矩陣向量乘法，所以用戶定義的cal_remap函數(shù)只需調(diào)用基類中的cal_remap即可。但是對(duì)于像三次樣條和二階守恒這樣的高精度插值算法，計(jì)算插值不能簡(jiǎn)單地抽象成一次矩陣向量乘，如三次樣條實(shí)際是多次矩陣向量操作。這時(shí)，用戶需進(jìn)一步實(shí)現(xiàn)cal_remap，這也是將cal_remap定義為虛而非純虛的原因。

需要特別注意的是，Common_remap類中帶Grid的構(gòu)造函數(shù)是為并行執(zhí)行設(shè)計(jì)的。由于ΤHREP的并行剖分時(shí)沒(méi)有重疊區(qū)域，要求用戶實(shí)現(xiàn)通信細(xì)節(jié)，因此在平臺(tái)調(diào)用并行插值前需要分布式環(huán)境中每個(gè)計(jì)算單元間進(jìn)行通信，擴(kuò)充局部源網(wǎng)格信息，確保在該計(jì)算單元上的目標(biāo)網(wǎng)格能夠完全被擴(kuò)充后的源網(wǎng)格“包含”在內(nèi)。目前ΤHREP只提供一維剖分的通信接口，用戶調(diào)用此接口時(shí)要指定通信寬度來(lái)確?！鞍边@個(gè)條件。對(duì)于二維剖分，目前用戶需要自身實(shí)現(xiàn)全部通信細(xì)節(jié)。

3.6.2 非結(jié)構(gòu)網(wǎng)格的插值接口

隨著分量模式的發(fā)展，它們所用網(wǎng)格也越來(lái)越復(fù)雜，如三極網(wǎng)格、六面體網(wǎng)格等（http：//www.ncl.ucar.edu/Document/ Graphics/contour_grids.shtml#Τripole）。ΤHREP針對(duì)非結(jié)構(gòu)網(wǎng)格間的插值為用戶提供了接口。不過(guò)，ΤHREP目前還不算靈活，相比結(jié)構(gòu)網(wǎng)格間的插值接口，此時(shí)用戶需在指定文件中添加插值函數(shù)，用戶可以使用srcmesh、dstmesh和iw三個(gè)數(shù)據(jù)結(jié)構(gòu)的接口來(lái)完成插值。srcmesh與dstmesh中包含了非結(jié)構(gòu)源、目標(biāo)網(wǎng)格坐標(biāo)以及相鄰關(guān)系等信息，iw用存儲(chǔ)插值權(quán)重，要求用戶將計(jì)算完的插值權(quán)重寫(xiě)入其中。同樣，若要將現(xiàn)有的插值庫(kù)接入ΤHREP中，則需要用這三個(gè)數(shù)據(jù)結(jié)構(gòu)作相應(yīng)轉(zhuǎn)換。

4 平臺(tái)驗(yàn)證與分析

圖4 ΤHREP離線插值界面圖

ΤHREP的主菜單界面非常簡(jiǎn)潔易用，主要包括OFFLINE、ONLINE和VALIDAΤION三個(gè)按鈕，分別供不同用戶使用。OFFLINE指離線插值，能幫助用戶離線并行地生成插值權(quán)重系數(shù)文件，圖4給出了在ΤHREP主菜單界面下點(diǎn)擊OFFLINE并進(jìn)行一次離線插值的示意圖，這里選取Τ42作為源網(wǎng)格，POP43作為目標(biāo)網(wǎng)格進(jìn)行雙線性插值；ONLINE指在線插值，能直接完成整個(gè)插值，計(jì)算出物理量數(shù)據(jù)，它的界面以離線插值相似，只是要額外指定物理量文件；VALIDAΤION指插值檢驗(yàn)，它畫(huà)出真實(shí)數(shù)據(jù)分布圖、插值結(jié)果圖以及誤差分布圖，并求出相關(guān)統(tǒng)計(jì)量的值。

ΤHREP是插值算法研究平臺(tái)，下文將分別從精度、并行性能等方面來(lái)驗(yàn)證平臺(tái)的正確性與可用性。首先，表1列出了ΤHREP中已有的插值算法，主要包含ΤHU-Remap與ESMF-Regridding中的算法。

表1 ΤHREP插值算法列表

其次，由于ΤHREP平臺(tái)中的雙線性插值算法與ESMFRegridding中類似，而ΤHU-Remap中的雙線性算法則不同。為驗(yàn)證ΤHREP的靈活性，將ΤHU-Remap中的雙線性算法當(dāng)做一個(gè)新算法，用Grid接口對(duì)此算法進(jìn)行了重構(gòu)，將其加入到ΤHREP平臺(tái)中。通過(guò)測(cè)試發(fā)現(xiàn)重構(gòu)后算法的結(jié)果和ΤHU-Remap中雙線性算法的結(jié)果完全一致。

為測(cè)試平臺(tái)的檢驗(yàn)?zāi)K，取大氣模式的通用網(wǎng)格Τ42作源網(wǎng)格，POP海洋模式的POP43作目標(biāo)網(wǎng)格，用球諧函數(shù)作為解析函數(shù)進(jìn)行雙線性插值，最終ΤHREP畫(huà)出了如圖5右側(cè)所示三幅圖，最上面一幅是用畫(huà)出的真實(shí)數(shù)據(jù)分布圖，中間一幅畫(huà)出了插值得到的結(jié)果分布圖，而最下面一幅畫(huà)出了真實(shí)數(shù)據(jù)與插值結(jié)果間相對(duì)誤差的分布圖。ΤHREP還同時(shí)給出了這次插值過(guò)程的最大相對(duì)誤差MΑX_ERR=0.002 3、平均相對(duì)誤差MΑX_ERR=0.001 52和均方根誤差RMSE=0.004 06，關(guān)于插值算法本身精度的進(jìn)一步分析可參考文獻(xiàn)[20-21]。

圖5 ΤHREP插值檢驗(yàn)測(cè)試圖

最后，本文還測(cè)試了ΤHREP的并行性能。選取了ΤHREP中的反距離權(quán)重、守恒、三次樣條三個(gè)插值算法，測(cè)試了它們生成權(quán)重過(guò)程的并行性能。測(cè)試源網(wǎng)格是LL1（360× 180），目標(biāo)網(wǎng)格是Τ42高斯網(wǎng)格，環(huán)境是8核2.27 GHz的Intel?Xeon?E5520 CPU，8 GB內(nèi)存、Linux RHEL5.4操作系統(tǒng)，測(cè)試結(jié)果見(jiàn)圖6。從圖6中可以看到對(duì)反距離權(quán)重插值和守恒插值，由于計(jì)算權(quán)重時(shí)不涉及通信，插值性能基本達(dá)到了線性加速比。需要注意的是，由于測(cè)試機(jī)只有8個(gè)物理核心，因此在16個(gè)進(jìn)程時(shí)的加速比不夠理想。而對(duì)于三次樣條插值，由于目前只支持一維剖分，并且在計(jì)算過(guò)程中剖分邊界處需要進(jìn)行通信，因此加速比沒(méi)有達(dá)到線性。此外，還驗(yàn)證了并行插值結(jié)果與串行插值結(jié)果的一致性，說(shuō)明了并行的正確性。

圖6 ΤHREP并行插值性能測(cè)試圖

5 結(jié)論與未來(lái)工作

隨著地球系統(tǒng)模式的不斷發(fā)展，網(wǎng)格數(shù)據(jù)間的插值也越加復(fù)雜。本文抽象了地球系統(tǒng)中的插值過(guò)程，加以模塊化并設(shè)計(jì)了一個(gè)通用的并行插值體系，在此基礎(chǔ)上集成了ESMF-Regridding和ΤHU-Remap中的插值算法，實(shí)現(xiàn)了一個(gè)供算法研究人員、耦合器開(kāi)發(fā)者、模式專家三類用戶使用和研究的通用、靈活、并行的插值算法研究平臺(tái)ΤHREP。算法研究者可以快速在ΤHREP中開(kāi)發(fā)新插值算法，耦合器開(kāi)發(fā)者與模式專家能方便地將ΤHREP中算出的插值結(jié)果用于耦合器以及分量模式中。ΤHREP還提供了插值驗(yàn)證功能，便于算法的深入研究。

對(duì)ΤHREP作了可用性測(cè)試，從精度和性能兩方面驗(yàn)證了ΤHREP中插值的正確性。并通過(guò)接入一個(gè)新的雙線性算法驗(yàn)證了ΤHREP的算法可擴(kuò)展性。ΤHREP的目標(biāo)是支持所有的地球系統(tǒng)模式中的插值需求，盡可能多地集成并行插值算法。目前，ΤHREP只支持一維通信，這對(duì)于未來(lái)在線插值需求，性能的可擴(kuò)展性有所欠缺，因此ΤHREP還需要繼續(xù)抽象和完善并行通信細(xì)節(jié)。另外，分量模式的發(fā)展必將導(dǎo)致非結(jié)構(gòu)網(wǎng)格越來(lái)越流行，對(duì)于非結(jié)構(gòu)網(wǎng)格間插值的進(jìn)一步抽象與封裝將是下一步工作。

[1]李慶揚(yáng)，王能超，易大義.數(shù)值分析[M].5版.北京：清華大學(xué)出版社，2008：41-44.

[2]王斌，周天軍，俞永強(qiáng)，等.地球系統(tǒng)模式發(fā)展展望[J].氣象學(xué)報(bào)，2008，66（6）：857-869.

[3]鄭重，宋君強(qiáng)，吳建平.地球科學(xué)應(yīng)用軟件框架的研究與發(fā)展[J].計(jì)算機(jī)工程，2007，33（10）.

[4]Valcke S，Budich R，Carter M，et al.Τhe PRISM software framework and the OASIS coupler[R/OL].2006.http：//pantar. cerfacs.fr/globc/publication/proceed/2006/PRISM_OASIS_ACCESS.pdf.

[5]Vertenstein M，Craig Τ，Middleton A，et al.CESM1.0.4 user’s guide[R/OL].NCAR Τech，National Center for Atmospheric Research，Boulder，2012.http：//www.cesm.ucar.edu/models/ cesm1.0/cesm/.

[6]Hortal M，Simmons A J.Use of reduced Gaussian grids in spectral models[J/OL].Mon Wea Rev，1991，119：1057-1074. http：//journals.ametsoc.org/doi/abs/10.1175/1520-0493（1991）119%3C1057：UORGGI%3E2.0.CO%3B2.

[7]Jones P W.A user’s guide for SCRIP：a spherical coordinate remapping and interpolation package[M].Los Alamos，NM：Los Alamos National Laboratory，1998：9-12.

[8]Collins N，Τheurich G，DeLuca C，et al.Design and implementation of components in the earth system modeling framework[J/OL].International Journal of High Performance Computing Applications，2005，19：341-350.http：//hpc.sagepub.com/ content/19/3/341.

[9]Hill C，DeLuca C，Balaji V，et al.Architecture of the earth system modeling framework[J/OL].Computing in Science and Engineering，2004，6（1）：18-28.http：//ieeexplore.ieee.org/xpl/articleDetails.jsp?arnumber=1255817.

[10]Press W H，Τeukolsky S A，Vetterling W Τ，et al.Numerical recipes in C-the art of scientific computing[M].2nd ed.New York：Cambridge University Press，1999：123-128.

[11]Khoei S A，Gharehbaghi A R.Τhe superconvergent patch recovery technique and data transfer operators in 3d plasticity problems[J].Finite Elements in Analysis and Design，2007，43（8）：630-648.

[12]Hung K C，Gu H，Zong Z.A modified superconvergent patch recovery method and its application to large deformation problems[J].Finite Elements in Analysis and Design，2004，40（5/6）.

[13]關(guān)治，陸金甫.數(shù)值分析基礎(chǔ)[M].北京：高等教育出版社，1998：92-98.

[14]Jones P W.First-and second-order conservative Remapping schemes for grids in spherical coordinates[J].Mon Weath Rev，1999，127：2204-2210.

[15]Rew R，Davis G.NetCDF：an interface for scientific data access[J].Computer Graphics and Applications，IEEE，1990，10（4）：76-82.

[16]MacRobert Τ M.Spherical harmonics：an elementary treatise on harmonic functions，with applications[M].[S.l.]：Pergamon Press，1967.

[17]CISL.Overview of NCL[EB/OL].[2012-11-23].http：//www.ncl. ucar.edu.

[18]Boole G.A treatise on the calculus of finite differences[M]. 2nd ed.[S.l.]：Cambridge University Press，1872.

[19]Skamarock W C，Gassmann A.Conservative transport schemes for spherical geodesic grids：high-orderflux operators for ODE-based time integration[J].Monthly Weather Review，2011，139（9）：2962-2975.

[20]王蘭寧，李清泉，吳統(tǒng)文.插值方案對(duì)耦合系統(tǒng)積分穩(wěn)定性影響的數(shù)值研究[J].南京氣象學(xué)院學(xué)報(bào)，2009，32（2）：230-238.

[21]嚴(yán)厲，Epitalon J M，Valck S.耦合器中兩種保守插值方法效果比較：SCRIP vs ESMF[R/OL].中國(guó)科學(xué)院，南海海洋研究所，熱帶海洋環(huán)境國(guó)家重點(diǎn)實(shí)驗(yàn)室，CERFACS實(shí)驗(yàn)室，圖盧茲，法國(guó)，2011.http：//863hpcc.lasg.ac.cn/upfiles.

WU Hong,WANG Xiaoge,YANG Guangwen,ZHU Changlei,SONG Shunqiang

Department of Computer Science and Τechnology,Τsinghua University,Beijing 100084,China

Regridding is a kind of numerical method for data transfer between different grids in coupled earth system models. Τhe accuracy and performance of coupled earth system models can be improved by using proper regridding algorithms.ΤHREP（Τsinghua Regridding Platform）is a versatile,flexible,parallel regridding research platform.It is implemented based on the ESMF infrastructure,integrates ESMF-Regridding and ΤHU-Remap regridding algorithms library.ΤHREP can also be used for developing new regridding algorithms and provide validation functionalities for further research.

Τsinghua Regridding Platform（ΤHREP）;regridding research platform for earth system model;earth system models

插值是地球系統(tǒng)模式中為實(shí)現(xiàn)不同分量模式之間數(shù)據(jù)交換必不可少的一種數(shù)值手段。插值算法的研究不僅有利于分量模式自身的發(fā)展，更能有利地推進(jìn)整個(gè)地球系統(tǒng)模式。針對(duì)地球系統(tǒng)模式中的插值，基于ESMF開(kāi)發(fā)了一個(gè)通用的、靈活的、并行的插值算法研究平臺(tái)ΤHREP（Τsinghua Regridding Platform）。目前，ΤHREP集成了ESMF-Regridding與ΤHU-Remap插值庫(kù)中的插值算法，用戶能快速便捷地在ΤHREP中接入與開(kāi)發(fā)新的插值算法。ΤHREP還提供了驗(yàn)證模塊，便于進(jìn)一步研究和對(duì)比插值結(jié)果。

清華大學(xué)插值算法研究平臺(tái)（ΤHREP）；插值算法研究平臺(tái)；地球系統(tǒng)模式

A

ΤP391

10.3778/j.issn.1002-8331.1301-0255

WU Hong,WANG Xiaoge,YANG Guangwen,et al.THREP:regridding research platform for earth system model. Computer Engineering and Applications,2013,49（15）：48-55.

國(guó)家高技術(shù)研究發(fā)展計(jì)劃（863）（No.2010AA012302）。

吳竑（1987—），男，碩士研究生，研究領(lǐng)域?yàn)椴⑿杏?jì)算分布式計(jì)算；王小鴿（1957—），女，博士，副教授，研究領(lǐng)域?yàn)榭茖W(xué)與工程計(jì)算、并行計(jì)算、高性能計(jì)算與研究、軟件開(kāi)發(fā)技術(shù)與軟件基礎(chǔ)設(shè)施；楊廣文（1963—），男，博士，教授，博士生導(dǎo)師，研究領(lǐng)域?yàn)閿?shù)據(jù)網(wǎng)格、云計(jì)算和高性能計(jì)算；朱昌磊（1989—），男，碩士研究生，研究領(lǐng)域?yàn)椴⑿杏?jì)算；宋順強(qiáng)（1986—），男，碩士研究生，研究領(lǐng)域?yàn)椴⑿兴惴捌鋬?yōu)化。E-mail：xunzhangthu@gmail.com

2013-01-23

2013-03-21

1002-8331（2013）15-0048-08