非線性入滲模型求參中參數初始值確定方法研究

沈洪政，王仰仁，韓娜娜

(天津農學院水利工程學院，天津 300384)

0 引 言

參數擬合初始值的確定是非線性模型參數優化的基本問題。對于非線性模型參數擬合人們提出了多種方法[1-3]。李春友[4]等提出用優化方法計算Van Genuchten模型參數，擬合初始值參考有關文獻選取。由于所選初值有時對擬合值有較大的影響,因此在這種情況下就需要多次調初值才可獲得較為滿意的擬合效果。錢天偉[5]等用單純形調優法和單純形加速法作為優化手段，結合計算機程序進行參數優化求解，但其初始值是任意給定的。以上均是從參數優化求解方法角度進行研究，對于優化初始值如何選取，沒有給出普遍適用的方法。

在求解非線性模型問題時，通常利用最小二乘法進行參數優化擬合。對于有明確物理意義的參數，可以根據其物理意義進行初始值的確定。對于沒有明確物理意義的模型參數，通常采用均勻網格法[6]和外推內插法[7]等進行初始值的選取。對于以上方法，由于需要進行反復尋優，選取的初始點較多，過程復雜，使用不方便。因此，本文以擬合精度較高的非線性土壤入滲Green Ampt模型和Horton公式為例，提出采用泰勒級數法確定初始值，避免了初始值選取的煩瑣過程，提高了參數擬合的效率和精度。

1 方法介紹

1.1 目標函數

(1)

1.2 入滲模型介紹

1.2.1 Green-Ampt模型

Green-Ampt模型是Green W H和Ampt G A根據毛管理論提出的近似積水模型[8,9]。在土層均勻、初始土壤干燥的薄層積水入滲條件下，利用達西定律，經過化簡與推導，可得到Green-Ampt入滲公式[10]：

(3)

式中：i(t)為入滲率，mm/min；ic為穩定入滲率，mm/min；I(t)為累積入滲量，mm；t為入滲時間，min；b為參數，其中ic和b為待定參數。

1.2.2 霍頓(Horton)公式

霍頓(1940年)提出的入滲公式為：

i(t)=ic+(i0-ic)e-βt

(4)

(5)

式中：i0為初始入滲率，cm/min；β為衰減參數，其余符號意義同前，其中i0、ic、β為待定參數。

Green-Ampt模型和Horton公式具有較高的擬合精度，因而被得到廣泛的使用。由公式的形式易知，兩者均為非線性模型，其參數的擬合屬于非線性擬合問題。在利用最小二乘法(OLS)進行參數優化擬合時，兩模型中的參數i0和ic的初始值可根據其物理意義進行確定，而參數b及β由于沒有明確的物理意義，用外推內插法等方法確定時過程煩瑣。針對此問題，本研究提出利用泰勒級數法對公式進行展開，結合田間實測入滲數據，獲得未知參數b和β的初始值，再利用OLS進行擬合，即可得到最優入滲模型參數。

1.3 初始值的確定

1.3.1 參數i0和ic初始值的確定

參數i0和ic分別表示初始入滲率和穩滲率。初始入滲率可選取入滲試驗測定的最大入滲率值作為參數i0的優化初始值，穩滲率可選取入滲試驗測定的最小入滲率值作為參數ic的優化初始值。

1.3.2 參數b初始值的確定

參數b初始值可采用泰勒級數法確定。為了便于泰勒公式的展開，對于公式(3)將其改寫為以t為因變量的顯式函數：

(6)

(7)

取式(7)中前3項，代入式(6)中得：

(8)

1.3.3 參數β初始值的確定

(9)

取式(9)中前3項，代入式(5)中得：

(10)

由式(8)和式(10)知，通過泰勒級數法展開后，非線性模型轉化為線性模型，結合i0和ic的實際物理意義及田間實測入滲數據即可確定入滲模型參數的初始值，進而可以利用最小二乘法進行參數擬合?，F以天津農學院農田水循環試驗基地測得的田間入滲數據為例，驗證方法的合理性。

2 實例計算

2.1 試驗區概況

天津農學院農田水循環試驗基地位于天津市楊柳青鎮大柳灘村(116°57′E，39°08′N，海拔5.49 m)，全年平均氣溫11.6 ℃。試驗田總面積1 hm2，地下水位變幅在2.06～3.70 m之間，試驗區土壤理化性質見表1，土壤水分特征曲線參數[11]見表2和圖1。

2.2 測試方法

試驗在2011年分2次測試，分別在灌溉冬水的麥田、不灌冬水的麥田和棉花茬地進行(試驗區簡圖見圖2)。測試時間為1月6日(氣溫為-5.6 ℃，相應麥田測試點次編號分別記為1和2)和3月28日(氣溫為11.2 ℃，麥田測試點次編號分別記為4和5)；棉花試驗田只在3月28日測試，其測試點次編號記為3。

表1 試驗區土壤理化性質

表2 試驗區土壤水分特征曲線參數

圖1 試驗區土壤水分特征曲線圖

試驗采用雙環入滲儀測定土壤水分入滲過程[12]。內環橫截面積1 000 cm2，內環直徑35.7 cm，外環直徑60 cm，內外環高度均為30 cm。試驗前，將滲透環埋于試驗田，埋入土壤的深度為10 cm。內環加水保持5 cm水頭，采用手動方式向環內加水，利用量筒計量，使得內外環水位齊平。試驗開始后記錄每次灌入水量和相應的時間，其內環灌入水量即為滲入水總量，測試直到入滲率基本穩定為止，記錄內環水量全部滲完為止的時間。

圖2 試驗區簡圖(單位：m)

2.3 結果分析

2.3.1 b值與β值計算結果分析

根據初始值確定方法的介紹，以天津農學院水循環試驗基地在2011年1月6日和2011年3月28日測得的入滲數據為基礎，經過計算，給出了各測試點b和β平均值的計算結果(表3)以及各測試點b和β值隨時間變化的曲線關系圖(圖3、圖4)。

表3 各測試點b與β平均值計算結果表

由表3知，通過平均值法計算得到的Green-Ampt模型中參數b的值在各測試點之間變化相對較大，Horton公式中參數β的值在個測試點之間變化較小。從Green-Ampt模型和Horton公式的實際意義分析來看，當入滲時間t→∞時，入滲速率i趨于穩定入滲率ic，當入滲時間t不斷增加時，入滲速率i也趨于穩定入滲率ic。為滿足此變化規律，要求b值和β值均為正值，這與計算的結果相符，表明用取平均值的方法計算b值和β值是可行的。

圖3 各測試點b值隨時間變化曲線

圖4 各測試點β值隨時間變化曲線

由圖3和圖4知，在入滲試驗初期，由于初始表層土壤含水率及土壤吸力大小等變化較大，入滲速率變化較大，導致b值和β值變化較大。隨著入滲試驗的進行，入滲速率逐漸穩定，各測試點b值和β值變化也趨于平穩。β值基本上在橫軸附近上下變化，b值基本上在79上下變化。說明，采用泰勒級數法確定的b值和β值具有較好的穩定性，且β值較b值更穩定。

2.3.2 Green-Ampt模型及Horton公式參數擬合結果

以b和β求出的值作為初始值，以累積入滲量I(t)的實測值與模擬值誤差平方和最小為目標，利用Excel規劃求解工具，分別對Green-Ampt模型及Horton公式的參數進行擬合。

表4 Horton公式參數擬合結果表

表5 Green-Ampt公式參數擬合結果表

由表4和表5知：

(1)采用泰勒級數法確定初始值來描述入滲變化具有較高精度，相關系數R2均達到0.85以上，最大值達到了0.998 6，具有較好的相關性。但由于初始入滲率i0和穩滲率ic只是簡單選取測試數據的最大值和最小值，同時β和b值均采用測試數據的平均值計算得出，因此，需對參數初始值進行擬合。由擬合后的結果知，相關系數R2均達到了0.99以上，擬合精度得到顯著提高。

(2)對比Horton公式和Green-Ampt模型參數擬合結果知，Horton公式擬合參數均為正值，Green-Ampt模型擬合結果中穩定入滲率 出現了負值，與穩定入滲率 的實際物理意義不符。從均方差計算結果來看，Horton公式計算的各測試點之間均方差的值相差不大，均在1.63上下變化。Green-Ampt模型計算的各測試點之間均方差的值相差較大，最大值達到了7.31，最小值為0.76，變化較大。因此，Horton公式的參數比Green-Ampt模型的參數更加穩定，公式更具有適用性。

圖5 Green-Ampt模型各測試點累積入滲量實測值與模擬值隨時間變化曲線

圖6 Horton公式各測試點累積入滲量實測值與模擬值隨時間變化曲線

圖5和圖6給出了各測試點Green-Ampt模型和Horton公式累積入滲量模擬值與對應實測值隨時間變化的曲線。由圖可見，各測試點Green-Ampt模型和Horton公式累積入滲量的模擬值擬合精度均較高。說明，采用泰勒級數法確定的初始值進行參數擬合具有較高的準確性。另一方面，由于測試時溫度的不同，灌冬水的麥田測試點1和測試點5，其累積入滲量之間變化較大，不灌冬水的麥田測試點2和測試點4也有相同變化規律。而對于測試溫度相同的麥田測試點1、2(-5.6 ℃)和測試點4、5(11.2 ℃)而言，其累積入滲量實測值與模擬值隨時間變化趨勢相近，表明，溫度對土壤水分入滲過程有很大影響。

2.3.3 溫度對入滲結果的影響

通過圖5和圖6的分析結果可知，溫度對入滲結果有很大影響，溫度越高，入滲速率越大，反之越小。茲以參數穩定性較好的Horton公式為例，考慮溫度對入滲結果的影響，在式(4)和式(5)中引入溫度修正系數，

i(t)=KT[ic+(i0-ic)e-βt]

(12)

其中：

KT=ea(T-15)

(13)

式中：KT為溫度修正系數；T為氣溫，℃；a為參數，15為參考溫度，℃。對式(11)，將各測試點(1和5、2和4)在同一時間的實測入滲量數值之比與不同溫度(-5.6 ℃，11.2 ℃)下的溫度修正系數KT值之比建立等式關系[式(14)]，由此便可求得參數a的值(見表6)，由測試點1和5求得的a值為0.044 1，由測試點2和4求得a的值為0.050 7。參數a的值確定以后，再將不同溫度(-5.6 ℃，11.2 ℃)值代入式(11)中，即可求得溫度修正系數值。

(14)

式中：T1、T2分別為不同灌水地塊的氣溫；I1(t)、I2(t)分別為不同測試點的實測入滲量值。圖7和圖8分別表示不考慮溫度修正的入滲曲線及考慮溫度修正的入滲曲線。由圖7和圖8知，引入溫度修正系數以后，不同溫度的入滲曲線很好的重合到一起(參考溫度條件下)，顯著地改善了入滲曲線擬合精度。

表6 a值計算結果表

圖7 不考慮溫度修正的入滲曲線

圖8 考慮溫度修正的入滲曲線

3 結 語

本文以非線性土壤入滲Green-Ampt模型和Horton公式為例，提出采用泰勒級數法確定其模型參數初始值，同時考慮了溫度對入滲結果的影響。結果表明，在初始值擬合基礎上，相關系數均達到0.85以上，經最小二乘法擬合后，實測值與模擬值相關系數提高到0.99以上?？紤]溫度修正系數以后，不同溫度的入滲曲線很好的重合到一起(參考溫度條件下)，顯著地改善了入滲曲線擬合精度。

泰勒級數法具有更好的理論基礎，避免了非線性規劃求參中需要設置多個初始值的不足，可顯著提高計算效率，而且可以提高參數尋優的準確度，適用于無法線性化的非線性模型參數初始值的確定。但是，在擬合計算過程中Green-Ampt模型參數 值出現了負值情況，與其實際物理意義相矛盾，尚需進一步研究。

□

[1] 朱良君，張光輝，任宗萍. 4種土壤入滲測定方法的比較[J]. 水土保持學報，2012，32(6):163-168.

[2] Ravindra V. Kale & Bhabagrahi Sahoo. Green-Ampt Infiltration Models for Varied Field Conditions: A Revisit[J]. Water Resource Manage，2011,25:3 505-3 536.

[3] 張新民，王根緒，胡想全，等. 優化技術在土壤入滲參數計算中的應用[J]. 灌溉排水學報，2004，(6): 52-54.

[4] 李春友，任 理，李保國. 利用優化方法求算Van Genuchten方程參數[J]. 水科學進展，2001，(12): 473-476.

[5] 錢天偉，陳繁榮，杜曉麗，等. 一種推求van Genuchten方程參數的高性能優化方法[J]. 土壤學報，2004， (11): 973-975.

[6] 尚松浩. 水資源系統分析方法及應用[M]. 北京：清華大學出版社，2012:179-180.

[7] 張 瑩. 運籌學基礎[M]. 北京：清華大學出版社，2004:135-137.

[8] Green W H，Ampt G A. Studies on soil physics:I. Flow of air and water through soils[J]. AgricSci，1911，(4):1-24.

[9] 李 毅，王全九，邵明安，等. Green-Ampt入滲模型及其應用[J]. 西北農林科技大學學報，2007，(2):225-230.

[10] 雷志棟，楊詩秀，謝森傳. 土壤水動力學[M]. 北京：清華大學出版社，1989:130-131.

[11] 車 政，王仰仁，王永紅，等. 農田土壤水分特征曲線參數擬合及其剖面變異特性研究[J]. 灌溉排水學報，2016:1-5.

[12] 周鑫洋，王仰仁，劉群昌. 四個入滲公式對農田潮土的適應性研究[J]. 灌溉排水學報，2015，(8):96-100.