趙玉杰，韓松俊，張寶忠，張曉龍，趙風(fēng)華

（1．流域水循環(huán)模擬與調(diào)控國家重點實驗室 中國水利水電科學(xué)研究院，北京 100038；2．中國科學(xué)院遺傳與發(fā)育生物學(xué)研究所農(nóng)業(yè)資源研究中心 中國科學(xué)院農(nóng)業(yè)水資源重點實驗室河北省節(jié)水農(nóng)業(yè)重點實驗室，河北 石家莊 050022；3．中國科學(xué)院地理科學(xué)與資源研究所 生態(tài)系統(tǒng)網(wǎng)絡(luò)觀測與模擬重點實驗室，北京 100101）

1 研究背景

華北平原（112°48′E—122°45′E，32°00′N—40°24′N）是我國最大的農(nóng)業(yè)生產(chǎn)區(qū)，冬小麥-夏玉米一年兩熟的種植模式約占華北平原農(nóng)田總面積的36%[1]，其對灌溉的依賴加劇了水資源緊缺狀況，確定其水分的蒸發(fā)消耗對灌溉用水管理至關(guān)重要。華北平原冬小麥在每年十月的中、下旬播種，次年六月初收獲，夏玉米在每年六月中、下旬播種，同年九月底至十月初收獲[2]，在作物收割后至下一種作物出苗前的一段時期（即間作期），地表覆蓋接近為裸地[3]。陸面特性的季節(jié)變化同時對大氣狀況造成影響，特別是間作期植被的減少抑制蒸發(fā)，引起溫度增加、濕度減小[4-5]，甚至引起區(qū)域降雨和環(huán)流的異常[6]。華北平原兩熟制農(nóng)田的陸面和大氣特性的顯著季節(jié)變化為利用Penman-Monteith公式等側(cè)重陸面過程的方法進行蒸發(fā)估算帶來了很多困難[7-8]，如無法反映陸面特性變化對區(qū)域氣候的反饋作用，需要較多陸面信息確定參數(shù)等[9]。相比之下，蒸發(fā)互補關(guān)系通過大氣蒸發(fā)能力間接反映陸面狀況的變化，相應(yīng)的蒸發(fā)估算公式不需要詳細的陸面信息[10]，在華北平原兩熟制農(nóng)田蒸發(fā)估算中具有一定優(yōu)勢。

互補關(guān)系最早由Bouchet[11]在1963年提出，現(xiàn)已從對稱的線性關(guān)系發(fā)展為描述陸面蒸發(fā)與大氣蒸發(fā)能力之間非線性相互作用規(guī)律的一種本構(gòu)關(guān)系[12-13]，目前主要存在三種方法：韓松俊等[14-15]的 S型公式，Brutsaert[16]的多項式公式以及Crago和Szilagyi等[17-18]在多項式公式基礎(chǔ)上提出的對自變量進行動態(tài)標(biāo)度的方法。前兩種方法通過參數(shù)考慮植被、土壤水分等陸面特性的影響[19-20]，而在基于動態(tài)標(biāo)度的互補關(guān)系公式中，標(biāo)度系數(shù)被認為間接反映了陸面干濕狀況的影響[21]。雖然Brutsaert的多項式公式及在其基礎(chǔ)上的動態(tài)標(biāo)度方法在不同生態(tài)系統(tǒng)取得了良好的效果，但在華北平原兩熟制農(nóng)田，間作期陸面特性變化劇烈，并不確定相關(guān)公式是否能夠僅利用常規(guī)氣象數(shù)據(jù)探知陸面特性變化對蒸散發(fā)的影響，因此本文利用四個通量站數(shù)據(jù)對此進行驗證分析。

2 研究數(shù)據(jù)和方法

2．1 研究站點和數(shù)據(jù)本研究收集了位于華北平原兩熟制農(nóng)田的四個通量站（禹城、欒城、館陶和位山）（表1）日或者半小時尺度的氣象數(shù)據(jù)（包括凈輻射、土壤熱通量、溫度、風(fēng)速、水汽壓/相對濕度、大氣壓強）和通量數(shù)據(jù)（潛熱通量和感熱通量）。禹城站和欒城站的數(shù)據(jù)來自中國通量觀測研究聯(lián)盟（http：//chinaflux．org/），館陶站的數(shù)據(jù)來自國家青藏高原科學(xué)數(shù)據(jù)中心（http：//data．tpdc．a(chǎn)c．cn），位山站的數(shù)據(jù)來自于清華大學(xué)位山生態(tài)水文觀測站。

表1 研究中使用的農(nóng)田通量站點基本信息

2．2 廣義互補公式與韓松俊等[14-15]的無量綱形式的非線性互補公式類似，Brutsaert[22]將傳統(tǒng)線性互補關(guān)系拓展到無量綱形式：

式中：E為實際蒸發(fā)量；Epo和Epa分別為潛在蒸發(fā)量和表觀潛在蒸發(fā)量。式（1）的函數(shù)形式為具有一個參數(shù)的四次多項式，但一般采用三次多項式（被稱為基本式）：

式中：x＝Epo/Epa，y＝E/Epa。Epo/Epa被認為在 0和 1之間，但分別利用 Priestley-Taylor（P-T）[23]公式和Penman[24]公式確定Epo和Epa時，在Epo一定的情況下，即使極端干旱Epa并不一定趨于無窮大，Epo/Epa并不能趨于0。為了滿足式（2）的邊界條件要求，Crago和 Szilagyi[17-18]先后對自變量 Epo/Epa乘以一個系數(shù)，強迫其在極端干旱的情況下趨近于0，得到動態(tài)標(biāo)度的自變量：

式中表觀潛在蒸發(fā)Epa通過Penman公式計算：

式中：Ta為氣溫；Δ為對應(yīng)的飽和水汽壓-空氣溫度斜率；γ為濕度計常數(shù)；Rn為凈輻射；G為土壤熱通量；es（Ta）和 ea分別為Ta下的飽和水氣壓和實際水汽壓；f（u）為風(fēng)速函數(shù)。為Epa的最大值，可通過假設(shè)參考高度水汽壓為0時的Penman公式計算[18]，即：

式中：Tdry為假設(shè)的區(qū)域蒸發(fā)為0時的干旱環(huán)境的氣溫，可根據(jù)實際氣溫與水汽壓確定[25]。在動態(tài)標(biāo)度方法中，潛在蒸發(fā)Epo通過溫度修正的P-T公式計算：

式中：Δ（Tw）為假設(shè)的濕潤環(huán)境氣溫Tw（可通過根據(jù)波文比確定的完全濕潤環(huán)境下的表面溫度近似確定[26]）對應(yīng)的飽和水汽壓-溫度曲線斜率；α為來自P-T系數(shù)的參數(shù)，最初采用1．26的默認值，目前被認為是反映陸面干濕狀態(tài)、植被等因素綜合影響的參數(shù)。

基于動態(tài)標(biāo)度的廣義互補公式在實踐中存在兩種不同的處理。Szilagyi等[18]直接將標(biāo)度后的自變量 X帶入 Brutsaert[22]的多項式公式：

Szilagyi等[18]提出了利用氣象數(shù)據(jù)提前確定參數(shù) α（約為1．15）的方法，得到了免參數(shù)校正的廣義互補蒸發(fā)估算方法[21-27]。

Crago和Qualls[17]則通過對多個通量站數(shù)據(jù)的分析認為線性函數(shù)更合理：

對于Epo、Crago和Qualls[17]認為不同站點的參數(shù)α具有差異，需要通過率定取得。

對于 Epa和中的風(fēng)速函數(shù)也存在兩種計算方法，Szilagyi[18]在研究中采用 Penman[24]在 1948年提出的風(fēng)速函數(shù)：

式中u2為地面上2 m高處的平均風(fēng)速。而Crago和Qualls[17]采用基于 Monin-Obukhov（M-O）相似理論的風(fēng)速函數(shù)：

式中：k＝0．41是 Karman常數(shù)；Rd＝0．287（KJ/（kg·k））為空氣的摩爾氣體常量；u為高度 z1處的風(fēng)速，m/s；z1為風(fēng)速的測量高度，m；z2為溫度和濕度的測量高度，m；d0為零平面位移高度，m；zov為控制動量傳輸?shù)拇植陂L度，m；zo為控制水熱傳輸?shù)牟诼书L度，m；d0＝0．67hc，zov＝0．123hc，zo＝0．1zov，hc為植被冠層高度。需要指出的是，實際的風(fēng)速函數(shù)隨著粗糙度和近地面空氣動力學(xué)狀況的變化具有動態(tài)變化，但互補關(guān)系的Epa中的風(fēng)速函數(shù)可視為一種參考風(fēng)速函數(shù)[28]，可采用固定的hc確定，在本研究中hc取植被的平均高度1．5 m。

本研究分析對比Brutsaert的多項式公式的基本式、標(biāo)度后的多項式和標(biāo)度后的線性公式在華北平原兩熟制農(nóng)田的適用性，并分析兩種風(fēng)速函數(shù)對估算效果的影響。在計算中采用半月的時間步長，根據(jù)華北平原兩熟制農(nóng)田特點，將一年分為三個階段進行分析：從十月下半月到次年五月結(jié)束的小麥季，從七月下半月到九月結(jié)束的玉米季，包括六月、七月上半月和十月上半月的間作期。在數(shù)據(jù)的質(zhì)量控制中，我們剔除了明顯不合理的數(shù)據(jù)，包括可用能量（Rn-G）和實際蒸發(fā)E小于0的數(shù)據(jù)、E/Epen大于1．2的數(shù)據(jù)。研究中使用絕對偏差MAE、均方根誤差RSME和納什效率系數(shù)NSE作為評價互補函數(shù)公式對實際蒸發(fā)量模擬效果的指標(biāo)，其中納什效率系數(shù)NSE可表示為：

3 結(jié)果與討論

3．1 采用Penman風(fēng)速函數(shù)的不同公式模擬效果圖1中給出了采用Penman風(fēng)速函數(shù)并設(shè)定參數(shù)α為1．26時各站點蒸發(fā)比E/Epen與標(biāo)度前后自變量的散點圖。散點分布于式（2）的多項式曲線和式（8）的直線左右，具有一定季節(jié)特性，其中小麥季的散點分布偏左，而玉米季和間作期的數(shù)據(jù)點大部分位于線的右側(cè)，且玉米季散點多位于右上部。散點的季節(jié)性分布反映了四個兩熟制站點的大氣環(huán)境受陸面特性變化的影響。三種方法都能較好的模擬實際蒸發(fā)量（表2），其中Brutsaert的基本式模擬效果最好，MAE介于 0．25～0．60 mm/d之間，均值為 0．42 mm/d，NSE在 0．73～0．92之間，均值為 0．84；標(biāo)度后的線性公式模擬效果略差于基本式，MAE在0．25～0．63 mm/d之間，均值為 0．44 mm/d，NSE在 0．68～0．92之間，均值為0．82；而標(biāo)度后的多項式公式模擬效果最差，MAE在0．37～0．75 mm/d之間，均值為0．54 mm/d，NSE在 0．57～0．84之間，均值為 0．74。

圖1 采用Penman風(fēng)速函數(shù)并設(shè)定α＝1．26時四個站點標(biāo)度前后的廣義互補關(guān)系散點圖

表2 采用Penman風(fēng)速函數(shù)（式（9））并分別固定α和率定參數(shù)α?xí)r三種公式對實際蒸發(fā)的模擬效果

通過最小化MAE對參數(shù)α進行率定（表2），基本式的參數(shù)α在禹城站為1．35，其他站點在1．23～1．24之間，標(biāo)度后的線性公式的參數(shù)α在1．20～1．30之間，小于基本式的率定值，而標(biāo)度后的多項式公式參數(shù) α介于1．26～1．60之間，在位山站率定值1．60明顯大于另外兩種公式的率定值。除部分站點，三種方法率定參數(shù)α都在1．26附近，表明α＝1．26在華北平原兩熟制農(nóng)田具有較強的適用性。

3．2 采用M-O風(fēng)速函數(shù)的不同公式模擬效果圖2顯示，固定參數(shù) α＝1．26并采用 M-O風(fēng)速函數(shù)時，各站點散點分布仍具有季節(jié)性特征，但三個公式的模擬效果都差于使用Penman風(fēng)速函數(shù)（表3）?；臼降?MAE增加至 0．35～1．34 mm/d，均值增加了 0．34 mm/d，NSE減小至-0．18～0．84，均值減小了0．37；標(biāo)度后的線性公式 MAE增加至 0．25～0．78 mm/d，均值增加了 0．10 mm/d，NSE減小到 0．53～0．91，均值減小了 0．10；標(biāo)度后的多項式公式的 MAE增加至 0．55～1．50 mm/d，均值增加 0．41 mm/d，NSE減小至-0．41～0．61，均值減小0．53。采用M-O風(fēng)速函數(shù)時，基本式在禹城站的模擬效果、標(biāo)度后多項式公式在禹城站和位山站的模擬效果大幅降低，而標(biāo)度后的線性公式的模擬效果較好。對參數(shù)α進行率定后，采用M-O風(fēng)速函數(shù)的不同公式的模擬效果有一定提高（表3），但基本式在禹城站、標(biāo)度后的多項式在禹城站和位山站的模擬效果仍然較低。圖3顯示，固定參數(shù)α＝1．26與整體率定α估算的實際蒸發(fā)都不能抓住五月的實際蒸發(fā)峰值。禹城站率定的α＝1．49（估算的平均日蒸發(fā)量為1．94 mm/d），固定參數(shù) α＝1．26（估算的平均日蒸發(fā)量為 1．61 mm/d）低估 30．34%。當(dāng)采用α＝1．26時，基本式和標(biāo)度后多項式公式分別低估了54．15%和64．61%的實際蒸發(fā)。在欒城、館陶和位山站，整體率定的α值與1．26非常接近，而分階段率定的 a值在小麥季遠大于 1．26（分別為 1．70、1．65和1．86），而在間作期小于 1．26（分別為 1．09、1．19和 1．00）。固定參數(shù) α＝1．26引起對實際蒸發(fā)的估算在小麥季低估 29．28%～38．97%，在間作期高估6．29%～24．26%，無法準(zhǔn)確反映蒸發(fā)量的季節(jié)變化。相比之下，采用階段性參數(shù)時，公式能夠準(zhǔn)確反映實際蒸發(fā)量在小麥快速生長季的峰值以及間作期的谷值。

圖2 采用M-O風(fēng)速函數(shù)并設(shè)定α＝1．26時四個站點標(biāo)度前后的廣義互補關(guān)系散點圖

表3 采用M-O風(fēng)速函數(shù)（式（10））并分別固定α和率定參數(shù)α?xí)r三種公式對實際蒸發(fā)的模擬效果

圖3 采用不同參數(shù)時標(biāo)度后線性公式估算的實際蒸發(fā)與觀測值的對比

Penman公式（式（4））是對表觀潛在蒸發(fā)的近似估算，對風(fēng)速函數(shù)的使用沒有確定的結(jié)論。M-O風(fēng)速函數(shù)需要通過冠層高度確定動量粗糙度、零平面位移高度和標(biāo)量粗糙度，這些參數(shù)具有季節(jié)性變化特征，在研究中我們采用年均冠層高度，未考慮冠層高度的季節(jié)性變化，這可能是引起公式在禹城和位山站估算效果較差的原因。Penman風(fēng)速函數(shù)適用于小到中等尺度粗糙度的表面，不需要冠層高度作為參數(shù)。對于缺乏表面粗糙度詳細信息的情況下，采用Penman風(fēng)速函數(shù)能夠合理的估算實際蒸發(fā)，但是存在E＞Epen的不合理的現(xiàn)象。

3．3 參數(shù)α的季節(jié)變化年內(nèi)使用相同參數(shù)時三個階段數(shù)據(jù)散點的分布表明廣義互補函數(shù)公式需要使用不同參數(shù)來捕捉季節(jié)的變化。在研究中分別率定小麥季、玉米季和間作期不同階段的參數(shù)α，無論采用Penman初始風(fēng)速函數(shù)還是M-O風(fēng)速函數(shù)，線性公式 y＝X和多項式公式y(tǒng)＝（2-x）x2對散點的擬合都更好（見圖4和圖5），但與圖5相比，采用Penman風(fēng)速函數(shù)存在更多自變量和因變量大于1的不合理情況，因此，后續(xù)分階段率定參數(shù)時采用M-O風(fēng)速函數(shù)。

圖4 四個通量站采用Penman風(fēng)速函數(shù)分階段率定參數(shù)α的廣義互補關(guān)系散點圖

圖5 四個通量站采用M-O風(fēng)速函數(shù)分階段率定參數(shù)α的廣義互補關(guān)系散點圖

采用M-O風(fēng)速函數(shù)分階段擬合參數(shù)α后，三種公式都能較好的估算所有站點的實際蒸發(fā)，對全年實際蒸發(fā)的模擬效果均有較大提高（表4），除禹城站外各公式的MAE減小至0．15～0．50 mm/d，0．10～0．46 mm/d，NSE增加至 0．80～0．97之間，增加了 0．06～0．54。圖6顯示，采用階段性參數(shù)時，公式能夠抓住實際蒸發(fā)在小麥快速生長季的峰值以及間作期的谷值。

表4 采用M-O風(fēng)速函數(shù)（式（10））時各公式分階段率定的參數(shù)值及其模擬效果

圖6 分階段率定參數(shù)對三種公式估算的實際蒸發(fā)與觀測值的對比

三種公式中，標(biāo)度后線性公式的模擬效果最優(yōu)?；臼胶蜆?biāo)度后線性公式在蒸發(fā)比較小時低估，在蒸發(fā)比較大時高估。在禹城站，當(dāng)實際蒸發(fā)較小時（十月至次年四月、六月），基本式和標(biāo)度后多項式公式對實際蒸發(fā)的低估量分別為30%和38%，當(dāng)實際蒸發(fā)較大時，兩種公式對實際蒸發(fā)分別高估33%和36%。標(biāo)度后線性公式（8）能夠較好反映實際蒸發(fā)的季節(jié)變化，估算效果優(yōu)于基本式和標(biāo)度后多項式公式，MAE分別降低 0．02～0．26 mm/d和 0．12～0．40 mm/d，而 NSE分別增加了 0～0．16和 0．07～0．31。

互補原理最初期望利用常規(guī)氣象數(shù)據(jù)能夠準(zhǔn)確及時地捕捉到陸面狀況的變化[29]，這樣可以采用固定的α值。但僅大氣狀況的變化可能并不能完全反映陸面狀況對蒸散發(fā)過程的影響[13]，且受到大尺度氣候條件變化的影響，造成參數(shù)α并不固定，受陸面水分條件和植被狀況等因素的影響而發(fā)生變化[30]。楊漢波等[31]的研究也表明互補關(guān)系中的參數(shù)α具有顯著的季節(jié)性。在華北平原兩熟制農(nóng)田，小麥季、玉米季和間作期的水分條件和植被狀況具有顯著差異，因此參數(shù)α明顯不同，而分階段率定效果更好，這也表明三種公式都可通過參數(shù)的變化來適應(yīng)環(huán)境的變化。

4 結(jié)論

本文以華北平原四個冬小麥-夏玉米兩熟制農(nóng)田通量站為研究對象，分析了Brutsaert的廣義互補多項式公式的基本式以及基于動態(tài)標(biāo)度的線性和多項式公式的適用性。主要得出了以下結(jié)論：

（1）采用Penman風(fēng)速函數(shù)時，三種公式都能有效地模擬華北平原兩熟制農(nóng)田的實際蒸發(fā)量，公式中的參數(shù)可設(shè)定為α＝1．26，Brutsaert的基本式對實際蒸發(fā)的估算效果最好。

（2）標(biāo)度后的線性公式受風(fēng)速函數(shù)影響最小，且適用于所有兩熟制農(nóng)田站點，而基本式和標(biāo)度后多項式公式在采用M-O風(fēng)速函數(shù)不能有效的反映實際蒸發(fā)的季節(jié)變化。

（3）分小麥季、玉米季和間作期分別擬合參數(shù)α不僅能夠顯著提高三種廣義互補公式對實際蒸發(fā)量的估算精度，同時可以縮小不同公式對實際蒸發(fā)量模擬的差距。