董勤各，許 迪，章少輝，白美健，李益農(nóng)

?

考慮灌溉參數(shù)空間變異的區(qū)域畦灌模擬與驗證

董勤各1,2，許 迪3※，章少輝3，白美健3，李益農(nóng)3

（1. 西北農(nóng)林科技大學水土保持研究所，楊凌 712100；2. 中國科學院水利部水土保持研究所，楊凌 712100；3. 中國水利水電科學研究院流域水循環(huán)國家重點實驗室，北京 100038）

精確評估區(qū)域畦田灌水質(zhì)量有助于提高農(nóng)田灌水管理水平，而具有空間變異性的灌溉參數(shù)如何有效表征是影響區(qū)域畦田灌水質(zhì)量精確模擬評價的關(guān)鍵因素。為此，該研究的目的在于借助Monte-Carlo抽樣，建立考慮畦灌參數(shù)空間變異性的區(qū)域畦灌模擬方法。采用Monte-Carlo抽樣將具有空間變異性的區(qū)域灌溉參數(shù)（如入畦單寬流量、土壤砂粒含量、黏粒含量、土壤容重等）離散表征為若干個灌溉參數(shù)樣本，依次輸入田塊尺度畦灌地表水流-土壤水動力學耦合模型，以模擬評價區(qū)域畦灌過程?；?次區(qū)域畦灌試驗的實測數(shù)據(jù)和1個對比的確定性畦灌模擬方法，驗證建立的模型的模擬效果。結(jié)果表明，所建模擬方法與確定性模擬方法模擬計算的灌水效率和灌水均勻具一定差異，所建模型的模擬值與實測值間的灌溉定額和田間水利用系數(shù)相對誤差分別為9.95%～12.23%和8.39%～10.21%，而基于現(xiàn)有模型的相對誤差則分別為14.15%～16.78%和13.87%～15.88%，畦田平均土壤含水率實測值與所建模型模擬值的累積分布趨勢表現(xiàn)出良好的一致性。上述指標表明所建模擬方法有效縮小了區(qū)域灌溉參數(shù)空間均化處理所帶來的模擬誤差范圍較大等問題，為區(qū)域畦田灌溉優(yōu)化設(shè)計和管理提供了技術(shù)支撐。

灌溉；Monte-Carlo方法；土壤; 參數(shù)；空間變異；灌水質(zhì)量；數(shù)值模擬

0 引 言

畦田灌溉運行管理簡單、省時節(jié)力成本低，占據(jù)著中國農(nóng)業(yè)灌溉方式的主導地位[1]。精確模擬評價畦田灌水質(zhì)量，對提高農(nóng)業(yè)用水管理水平、促進灌溉農(nóng)業(yè)的可持續(xù)發(fā)展具有重要意義?；谖锢頇C制的畦灌模擬模型（如WinSRFR與B2B等）已在畦灌過程模擬與灌水質(zhì)量評價方面得到廣泛應用[2-3]。這些畦灌模擬模型多用于模擬單個田塊灌溉過程[4-7]，難以描述區(qū)域畦灌參數(shù)空間變異下的灌溉水流動態(tài)變化規(guī)律，致使畦灌模擬評價往往存在難以解釋的誤差。針對土壤特性、入畦單寬流量、田面糙率、田面微地形等灌溉參數(shù)區(qū)域空間變異性的處理方式，區(qū)域農(nóng)田水文過程模擬模型可分為3類：1）確定性模擬模型。該類模型主要將灌溉參數(shù)在區(qū)域尺度上均化處理，得到1組確定性灌溉參數(shù)值，輸入田塊尺度農(nóng)田水文模擬模型進行模擬計算[8-9]；2）分布式模擬模型。該類模型將模擬區(qū)域劃分為若干子區(qū)域，對子區(qū)域內(nèi)灌溉參數(shù)加權(quán)平均，并依次將所得灌溉參數(shù)值輸入給田塊尺度農(nóng)田水文模型，然后對上述模擬數(shù)據(jù)依據(jù)子區(qū)域面積進行加權(quán)計算，實現(xiàn)對區(qū)域農(nóng)田水文過程的有效模擬[10-12]；3）空間隨機模擬模型。該類模型借助參數(shù)隨機抽樣方法，充分表征灌溉參數(shù)的區(qū)域空間變異性，然后將抽取的灌溉參數(shù)值依次輸入給田塊尺度農(nóng)田水文模型，對模擬數(shù)據(jù)進行同系數(shù)加權(quán)計算，以模擬區(qū)域農(nóng)田水文過程[13]。確定性模擬與分布式模擬均難以精確描述灌溉參數(shù)空間變異性對區(qū)域灌溉過程及灌水質(zhì)量的影響，且分布式模擬易導致子區(qū)域灌溉參數(shù)值被過高或過低估計[8-12]，嚴重影響模擬精度；參數(shù)隨機模擬能有效避免上述問題，然而鮮見其用于區(qū)域畦灌過程模擬方面的報道。

Monte-Carlo抽樣是一種以概率統(tǒng)計理論為主要理論基礎(chǔ)的參數(shù)隨機模擬方法，能夠避開研究對象復雜內(nèi)部特性描述及其在系統(tǒng)行為上的困難，直接以系統(tǒng)運用過程模擬替代系統(tǒng)分析，便于解決那些由于計算過于復雜而難以得到解析解或者根本沒有解析解的問題[14-17]。因此，針對現(xiàn)有灌溉模型參數(shù)過于受到人為因素干擾，而未能很好表征區(qū)域畦灌參數(shù)的空間變異性，從而影響模擬精度等問題，本文以已有的畦灌地表水流-土壤水動力 學耦合模擬模型為工具[7-8]，借助Monte-Carlo抽樣，建立考慮畦灌參數(shù)空間變異性的區(qū)域畦灌模擬方法（stochastic parameter irrigation model，SRM）。基于在山東麻灣灌區(qū)開展的3次區(qū)域畦田灌溉試驗，對模型的模擬效果進行驗證，以期為區(qū)域畦田灌溉優(yōu)化設(shè)計和管理提供技術(shù)支撐。

1 考慮灌溉參數(shù)空間變異的區(qū)域畦灌模擬

1.1 模型假設(shè)

本文建立的考慮灌溉參數(shù)空間變異性的區(qū)域畦灌模擬方法主要包括2個部分，即畦灌地表水流-土壤水動力學耦合模擬與灌溉參數(shù)的Monte-Carlo抽樣。為簡化區(qū)域畦灌模擬過程，特假設(shè)如下：1）單個田塊內(nèi)土壤質(zhì)地均勻，考慮到區(qū)域內(nèi)農(nóng)田耕作與管理措施相同，坡度可視為一致；2）畦田長寬比較大，地表水流橫向運動過程遠遠小于沿畦長方向的縱向運動過程，一般將畦灌地表水流為沿畦長方向的一維運動過程；3）由于土壤沿畦長方向的水平入滲速率遠遠低于地表水流運動速度，導致土壤垂向入滲占主導地位，故只模擬土壤垂向入滲過程；4）只模擬田塊內(nèi)畦灌水流過程，不涉及渠系水流過程，且田塊間水流通量忽略不計。

1.2 畦灌地表水流-土壤水動力學耦合模擬

本文采用一維畦灌地表水流運動模型模擬畦灌地表水流運動過程、一維土壤水動力學模型模擬灌溉水流入滲過程及土壤水再分布過程。一維畦灌地表水流運動模型采用基于有限差分法、有限體積法與有限單元法建立的混合數(shù)值方法解算，一維土壤水動力學模型采用高精度有限體積法與有限差分格式解算，具體解算過程詳見文獻[7,18]。

1.2.1 控制方程

單個田塊的畦灌地表-土壤水流運動過程可采用由一維地表水流運動模型和一維土壤水動力學模型耦合構(gòu)建的田塊尺度畦灌數(shù)值模型描述[7-8]。在假設(shè)單個畦田的土壤質(zhì)地均勻分布基礎(chǔ)上，畦灌一維地表水流運動控制方程的矢量形式可表達如下

U+F=1+2+3（1）

式中U是因變量矢量；F是物理通量；1是畦面糙率；2和3分別是地形項和入滲項，分別用于描述畦田微地形狀況及土壤特性對畦灌水流運動的影響。

，，（3）

式中h為地表水深，m；為單寬流量，m3/(s·m)；為入滲率，可通過式（4）計算土壤水分單位時間變化量獲得，m/s；為重力加速度，m/s2；為水流推進距離，m；z為地面高程，m；為時間，s；為曼寧糙率系數(shù)，m1/6；為垂直平均流速，m/s。

畦灌時地表水流沿畦長方向運動速度遠大于地表水流入滲后形成的土壤水流水平運動速度，且地表水流沿畦長方向推進的同時還垂向入滲土壤，故畦灌時土壤水流沿畦長方向運移的濕潤鋒位置取決于地表水流的推進鋒位置，則土壤水平入滲過程可忽略，僅考慮垂向入滲過程[4,7]。同時，根據(jù)試驗觀測發(fā)現(xiàn)，單個田塊畦灌過程約為0.5～1.5 h，一次區(qū)域畦灌過程約為3 d?？梢姡喔冗^程持續(xù)時間短、灌溉水流入滲及再分布過程在此期間受冬小麥根系吸水過程影響小，可忽略不計。因此，可將畦灌土壤水運移簡化成不考慮作物吸水過程的一維非飽和土壤水動力學問題。

式中為土壤壓力水頭，入滲模擬時可用h代替，cm；為垂向距離，向下為正，cm；()為土壤非飽和導水率，cm/min；()為比水容量；為時間，min。

式（4）中的()和()可分別表示如下

（6）

（7）

式中K為土壤飽和導水率，cm/min；S為土壤有效飽和度[7]；、和均為土壤水分特征曲線參數(shù)，=1?1/；θ和θ分別為土壤殘余及飽和含水率，cm3/cm3；為壓力水頭，cm；()為土壤含水率，cm3/cm3，其與關(guān)系采用van Genuchten模型[4]（VG模型）表達如下

1.2.2 初始與邊界條件

地表水流初始條件為畦田內(nèi)地表水深和平均均為0，即h=0和=0。邊界條件為：當畦首處于入流狀態(tài)時，給定單寬流量=0（0為灌溉開始后的入畦單寬流量）；灌水停止后，畦首單寬流量=0，由于灌溉地表水流始終處于亞臨界流狀態(tài)，故畦首水深值由畦內(nèi)水深值插值獲得。此外，當畦尾處于封口狀態(tài)時，=0；若畦尾處于敞口狀態(tài)，則依據(jù)灌溉水流始終處于亞臨界流狀態(tài)的物理事實[19]，采用畦內(nèi)單寬流量和水深值插值獲得畦尾單寬流量和水深h值。

土壤水流初始條件為

式中h()為土壤初始負壓水頭分布。

邊界條件包括上邊界條件與下邊界條件。當畦灌強度（單位面積單位時間上畦田灌溉水深）大于土壤入滲能力時，產(chǎn)生地表積水或徑流，此時上邊界條件為

式中h為地表積水深度，其值為100h，cm；()為時刻的地表蒸發(fā)強度與畦灌強度之差，cm/min；h為地表最大積水深度，cm；h為土壤表層初始壓力水頭，cm；其他符號物理意義同上。

當?shù)叵滤惠^深時，下邊界可采用自由排水條件

式中為下邊界深度，cm。冬小麥根系主要集中在0～100 cm土層內(nèi)，故取值為100 cm。

1.3 考慮空間變異性的畦灌參數(shù)隨機模擬

1.3.1 區(qū)域畦灌參數(shù)的統(tǒng)計分析

在上述一維畦灌地表水流-土壤水動力學耦合模型中，將畦田長度、入畦單寬流量、土壤水力特性、畦田微地形狀況作為初始與邊界條件輸入模型。由于土壤水力特性與土壤特性（如土壤砂粒、黏粒含量和有機質(zhì)含量）呈非線性函數(shù)關(guān)系，故可由土壤砂粒、粉粒、黏粒和有機質(zhì)含量耦合表達。在評價區(qū)域畦田灌水質(zhì)量過程中，由于土壤特性、入畦單寬流量、田塊規(guī)格等均具有空間不確定性特征，諸多學者常采用正態(tài)分布或?qū)?shù)正態(tài)分布描述其分布規(guī)律[13,20-22]。區(qū)域尺度上，畦田長度、畦田寬度、入畦單寬流量、土壤砂粒、粉粒、黏粒和有機質(zhì)含量均呈現(xiàn)近似正態(tài)分布特征[20]。正態(tài)分布的概率分布函數(shù)為

式中為隨機變量；為隨機變量的均值；為隨機變量的標準差。利用Monte-Carlo隨機生成的任意1組灌溉參數(shù)樣本的值域范圍通常在[–∞，+∞]之間，但由于農(nóng)田生產(chǎn)活動和灌溉管理的要求，灌溉參數(shù)值域范圍是有限的，而且與參數(shù)的標準偏差值密切相關(guān)。為了真實反映灌溉參數(shù)的變化幅度范圍，需對隨機生成的參數(shù)值域范圍進行必要的限制。從實際意義和數(shù)學角度考慮，區(qū)域畦灌灌溉參數(shù)的值域范圍可限制在[–3，＋3]之間[23]，這是由于該范圍為置信水平為95%以上的置信區(qū)域。

畦田微地形狀況包括畦田微地形起伏幅度及起伏位置空間分布[22]。由于相同耕作措施下的區(qū)域畦田微地形起伏幅度具有相似性，故量化畦田微地形起伏幅度的田面相對高程標準差可近似相同。然而，相同標準差下的畦田微地形起伏位置空間分布并非唯一，其差異性對灌水質(zhì)量存在影響，需考慮多個畦田微地形起伏位置空間分布狀況表征區(qū)域尺度上的畦田微地形狀況[22]。用于表征畦田微地形狀況的田面相對高程呈現(xiàn)近似正態(tài)分布特征[12]。

1.3.2 Monte-Carlo模擬

Monte-Carlo模擬以概率統(tǒng)計理論為基礎(chǔ)[23]，可用來分析各種不確定性問題。它假定隨機變量的概率分布函數(shù)是已知的，通過隨機抽樣方法得到隨機變量的若干個離散樣本，通過數(shù)值模擬得到相應的輸出結(jié)果，通過對輸出結(jié)果的統(tǒng)計，得到均值、方差等估計量，并擬合其輸出結(jié)果的概率分布特征。畦田灌水質(zhì)量受畦田長度、入畦單寬流量、土壤特性、畦田微地形狀況等灌溉參數(shù)約束。

本文假定區(qū)域畦田內(nèi)的坡度、田面糙率與微地形標準差為常值，土壤特性、畦田長度、入畦單寬流量、畦田微地形狀況的空間變異性是影響灌水質(zhì)量優(yōu)劣的重要因素。具體模擬流程如下：1）根據(jù)擬合統(tǒng)計得到的畦田長度、入畦單寬流量概率分布函數(shù)進行隨機抽樣，獲得表征兩者空間變異性的區(qū)域離散樣本；2）根據(jù)擬合統(tǒng)計得到的土壤粉粒、黏粒及容重的概率分布函數(shù)，對上述變量進行隨機抽樣，并將隨機抽取的樣本輸入Rosetta人工神經(jīng)網(wǎng)絡(luò)模型[24-25]，轉(zhuǎn)換為模型所需的土壤水力特性參數(shù)（飽和含水率θ、殘余含水率θ、飽和導水率K、參數(shù)、）；3）通過畦田長度、寬度及田塊尺度畦灌模型地表單元格長度，計算田面相對高程數(shù)量，借助田面微地形隨機模擬模型[22]，獲得田面相對高程序列，并對其沿畦寬方向做均化處理與統(tǒng)計特征修正，作為一個畦田微地形狀況樣本輸入給田塊尺度畦灌數(shù)值模型；4）重復上述步驟，獲得由畦田長度、入畦單寬流量、土壤水力特性及畦田微地形狀況組成的樣本，依次輸入畦灌地表水流-土壤水動力學耦合模型進行模擬；5）基于上述模擬結(jié)果，計算灌水質(zhì)量評價指標，評價區(qū)域畦田灌水質(zhì)量。區(qū)域畦灌隨機模擬模型的構(gòu)建流程如圖1所示。

1.3.3 灌水質(zhì)量評價指標

本文采用灌水效率與灌水均勻度作為灌水質(zhì)量評價指標[26-27]。灌水效率通常是指1次灌水存儲在作物根區(qū)的平均水深和畦田平均灌水深度之比。灌水均勻度是反映灌水沿畦田長度方向的入滲分布均勻程度。灌水均勻度越大，田塊內(nèi)各點土壤儲存水分的均勻性就越好。

灌水效率E的計算公式如下[27]：

（14）

式中Z為灌溉后儲存在作物根區(qū)的平均水深，cm；Z為單個田塊平均灌水深度，cm；θ和θ分別為灌前和灌后計劃濕潤層內(nèi)平均土壤體積含水率，cm3/cm3；R為作物根區(qū)深度，依據(jù)實測情況取值為100 cm。

灌水均勻度C的計算公式如下[27]

式中C為灌水均勻度，%；Z為第個測點灌水深度，cm；為單個畦田沿畦長方向的測點數(shù)。

1.3.4 數(shù)據(jù)分析與處理

通過Monte-Carlo抽樣與畦灌地表水流-土壤水動力學耦合模型的模擬計算，針對抽樣得到的每組輸入值，統(tǒng)計分析模擬輸出數(shù)據(jù)，計算灌水質(zhì)量評價指標。由于一維畦灌地表水流-土壤水動力學耦合模型的輸出為單個田塊的畦灌過程數(shù)據(jù)，而輸入不確定性為具有空間隨機性的土壤特性、畦田長度、入畦單寬流量與畦田微地形等灌溉參數(shù)。因此，為表征區(qū)域尺度畦灌參數(shù)不確定性的影響，本文采用累積概率分布及經(jīng)典統(tǒng)計方法，驗證考慮畦灌參數(shù)隨機模擬的區(qū)域畦灌模擬方法。

2 區(qū)域畦灌試驗與模型參數(shù)選取及樣本量確定

2.1 區(qū)域畦灌試驗

2.1.1 試驗區(qū)與觀測指標

田間試驗在山東東營市麻灣灌區(qū)開展（118°22′E，37°12′N），試驗區(qū)面積76 hm2，包括272塊畦田。麻灣灌區(qū)地處黃河下游濱海地區(qū)，屬暖溫帶季風型大陸性氣候，多年平均降水量為537.4 mm。試驗區(qū)內(nèi)0～100 cm土層以黏壤土為主。主要種植作物為冬小麥和夏玉米。該灌區(qū)主要在冬小麥越冬期前（11月）、返青-拔節(jié)期（3月）和抽穗-灌漿期（5月）進行灌溉，補充作物所需水分。一般采用條畦灌溉方式，水源為黃河引水。由于冬小麥灌漿期植株過高過密，不利于精確觀測水流推進過程和畦面微地形狀況，因此，分別選在2012年3月（T1），2012年11月（T2），2013年3月（T3）開展3次冬小麥區(qū)域畦灌試驗，以驗證所建區(qū)域畦灌模擬方法的準確性。

畦灌前，采用工程測量繩測定272塊畦田的長度與寬度；利用自動對焦水準儀（蘇州一光儀器有限公司，±1.0 mm）觀測選取的100塊畦田的相對高程。畦灌時，觀測100塊畦田的水流推進過程和入畦流量，水流推進到畦尾后關(guān)閉進水口，便攜式自動流量計（上海安鈞電子科技有限公司，±1%）測定入畦流量。灌前1 d與灌后2 d，在隨機選取的100塊畦田的畦首、畦中與畦尾處，均沿地面垂向分4層（0～20 cm、>20～40 cm、>40～70 cm及>70～100 cm）采集土樣測定土壤含水率，其中采用畦田中部前3層土樣測定土壤顆粒組成、容重和有機質(zhì)含量。采用比重計法測定土壤顆粒組成，重鉻酸鉀法測定土壤有機質(zhì)含量，烘干法測定土壤容重和土壤含水率[28]。

2.1.2 區(qū)域畦灌參數(shù)統(tǒng)計特征

畦田灌水質(zhì)量主要受畦田長度、田面微地形狀況、土壤特性與入畦單寬流量等灌溉參數(shù)影響。根據(jù)研究區(qū)內(nèi)畦田分布與耕作管理特性，畦田長度表現(xiàn)出相對穩(wěn)定的特征，其概率分布特征如圖2a所示。擬合統(tǒng)計可以看出，研究區(qū)內(nèi)畦田長度基本呈現(xiàn)正態(tài)分布特征，其概率密度函數(shù)可采用正態(tài)分布函數(shù)來擬合，擬合分布模型如式（12）所示。畦田長度正態(tài)擬合的參數(shù)有2個：均值120 m，標準差33.6 m。研究區(qū)內(nèi)畦田寬度具有一定的空間隨機性（圖2b），呈現(xiàn)近似正態(tài)分布統(tǒng)計特征，其擬合分布模型如式（12）所示，對研究區(qū)內(nèi)畦田寬度擬合，均值9 m，標準差1.98 m。考慮到模型模擬時采用的是入畦單寬流量，基于觀測的入畦流量以及畦田寬度數(shù)據(jù)，計算每塊畦田的入畦單寬流量。表1給出入畦單寬流量的統(tǒng)計特征和K-S檢驗結(jié)果，可以看出，3次試驗的入畦單寬流量平均值的變化很小，但相應的標準差變化幅度較大；經(jīng)Kolomogorov-Smirnov方法（K-S）檢驗[29]，入畦單寬流量呈現(xiàn)明顯的正態(tài)分布統(tǒng)計特征。

表1 麻灣灌區(qū)試驗區(qū)土壤容重和有機質(zhì)及單寬流量的統(tǒng)計特征及K-S檢驗結(jié)果

注：N為正態(tài)分布，下同。

Note: N is normal distribution, the same as below.

研究區(qū)采用相似的農(nóng)田耕作管理方式，隨機選取的100塊畦田的田面相對高程的標準差比較一致。3次試驗的田面相對高程標準差的均值及其標準差分別為4.06、4.23、4.17 cm與0.12、0.15、0.14 cm，畦面坡度均近似為0.6‰，田面相對高程空間變異特征參數(shù)均值如表2所示，可用于模擬田面微地形。由于畦灌地表水流-土壤水動力學耦合模型為沿畦長方向的一維數(shù)值模擬，因此需要對通過田面微地形隨機模擬模型[22]生成的田面相對高程數(shù)據(jù)沿畦寬方向做均化處理。

表2 麻灣灌區(qū)田面相對高程空間變異特征參數(shù)

土壤特性主要包括土壤質(zhì)地、容重及有機質(zhì)含量。根據(jù)測定的土壤質(zhì)地、容重與有機質(zhì)含量數(shù)據(jù)，分別統(tǒng)計其概率分布特征。經(jīng)K-S檢驗，試驗區(qū)內(nèi)每層土壤的砂粒、粉粒、黏粒含量均呈現(xiàn)明顯的正態(tài)分布統(tǒng)計特征，分布模型如式（12）所示。3次試驗測定發(fā)現(xiàn)，土壤質(zhì)地統(tǒng)計特征值基本一致，故本文采用2012年3月試驗的數(shù)據(jù)統(tǒng)計指標（表3）隨機抽取土壤質(zhì)地樣本。3次試驗的土壤容重均呈現(xiàn)明顯的正態(tài)分布統(tǒng)計特征，主要參數(shù)如表1所示。3次試驗的土壤有機質(zhì)含量均值與標準差比較一致，均呈現(xiàn)正態(tài)分布特征（表1）。

表3 麻彎灌區(qū)土壤顆粒組成的統(tǒng)計特征及K-S檢驗結(jié)果

依據(jù)美國農(nóng)業(yè)干旱研究中心研究成果，畦灌過程中田面糙率系數(shù)一般可參照作物長勢獲得[30-31]，故通過參考其提供的基于冬小麥生育期長勢確定的田面糙率系數(shù)數(shù)據(jù)庫，將麻灣灌區(qū)3次試驗的田面糙率系數(shù)分別設(shè)為0.13、0.09和0.14。

2.2 畦灌參數(shù)選取及參數(shù)樣本量確定方法

模型的輸入?yún)?shù)包括畦長、田面微地形狀況、入畦單寬流量、灌水時間、地表糙率系數(shù)和由土壤特性參數(shù)推求的VG模型參數(shù)（θ、θ、K、、）。畦長取實測值均值作為模型輸入?yún)?shù)；入畦單寬流量通過Monte-Carlo抽樣獲得；土壤水力特性參數(shù)由Monte-Carlo抽取土壤砂粒、粉粒及容重樣本后經(jīng)Rosetta模型轉(zhuǎn)換獲得；田面微地形狀況由修正的田面微地形隨機模擬模型生成并修正后獲得；3次試驗的地表糙率系數(shù)，依據(jù)美國農(nóng)業(yè)干旱研究中心研究成果分別設(shè)為0.13、0.09和0.14；依據(jù)當?shù)剞r(nóng)田灌水實際，水流推進到畦田尾部時關(guān)閉入水口。

采用Monte-Carlo抽樣研究區(qū)域灌溉參數(shù)空間變異性時，需依據(jù)各參數(shù)的概率分布函數(shù)，隨機抽取與研究區(qū)灌溉參數(shù)分布特征相一致的參數(shù)樣本量。將抽取的隨機樣本作為畦灌地表水流-土壤水動力學耦合模型的輸入條件。當隨機抽取一定數(shù)量的隨機參數(shù)組合樣本后，需要解決的問題是，針對給定模型輸入?yún)?shù)特征值，運行田塊尺度畦灌耦合模型評價區(qū)域畦田灌水質(zhì)量時，模擬多少組隨機參數(shù)樣本組合才能較真實表征實際區(qū)域畦灌過程、合理評價區(qū)域畦田灌水質(zhì)量，這需要設(shè)定科學嚴謹?shù)呐袛鄻藴?。在其他條件相同的情況下，給定的隨機參數(shù)樣本組合對應唯一的灌水效率E和灌水均勻度C。因此，從隨機參數(shù)分布總體中抽取一定數(shù)量的樣本作為代表，利用畦灌耦合模型進行畦灌模擬試驗，借助式（16）、（17）與（18）統(tǒng)計分析若干參數(shù)樣本條件下得到的區(qū)域畦田灌水質(zhì)量評價指標，并得出可以反映隨機參數(shù)分布總體對灌水質(zhì)量影響的最小個體數(shù)量。當模擬次數(shù)達到一定數(shù)量后，畦田灌水質(zhì)量評價指標值均值M和標準差R將趨于穩(wěn)定。如果以連續(xù)6個數(shù)據(jù)點的M<0.5%和R<5%作為均值和標準差達到穩(wěn)定狀態(tài)的判別指標，則可確定灌水質(zhì)量評價指標值均值及相應標準差。M和R的計算公式分別如下

（17）

（18）

式中M為第+1個灌水質(zhì)量評價指標值均值與第個的相對誤差；R為+1個灌水質(zhì)量評價指標值的標準差與個的相對誤差；X+1為基于第+ 1個參數(shù)樣本模擬計算的灌水質(zhì)量評價指標值；為基于第1～個參數(shù)樣本模擬計算的灌水質(zhì)量評價指標值的均值；S為個灌水質(zhì)量評價指標值的標準差；為灌水質(zhì)量評價指標均值收斂時的模擬計算次數(shù)。

3 考慮灌溉參數(shù)空間變異的區(qū)域畦灌模型驗證

3.1 對比的灌溉模擬方法

采用確定性模擬方法作為對比的灌溉模擬方法（deterministic parameter irrigation model，F(xiàn)IM），該方法是將模型輸入?yún)?shù)在區(qū)域尺度上進行平均處理，其值輸入給模型以進行模擬，并輸出單一灌溉過程數(shù)據(jù)，僅能描述區(qū)域整體的靜態(tài)灌溉特征。對現(xiàn)有模擬方法的改進，主要體現(xiàn)在模型輸入?yún)?shù)值如何精確描述區(qū)域灌溉參數(shù)空間變異性的差異上。對比的確定性模擬方法的輸入?yún)?shù)值如下：基于3次區(qū)域畦灌試驗確定的畦田長度均為120 m；糙率系數(shù)分別為0.13、0.09、0.14；入畦單寬流量分別為3.12、3.24、3.16 L/(s×m)；田面相對高程標準差的均值分別為4.06、4.23、4.17 cm（對應的空間特征變異參數(shù)見表2）；VG模型參數(shù)θ分別為0.428、0.442、0.421 cm3/cm3，θ分別為0.071、0.073、0.070 cm3/cm3，K分別為0.013 3、0.018 0、0.011 5 cm/min，分別為0.005 3、0.005 1、0.005 4，分別為1.658、1.667、1.651；依據(jù)當?shù)剞r(nóng)田灌水實際，水流推進到畦田尾部時關(guān)閉入水口。模擬每個畦灌過程時，畦長方向和土壤垂向的空間離散步長分別取為2.0和0.01 m，土壤與地表水流時間離散步長均取值為1 s。

3.2 驗證指標

選取灌溉定額、根系層土壤含水率均值與田間水利用系數(shù)作為驗證區(qū)域畦灌模型模擬精度的評價指標。灌溉定額是指進入到試驗區(qū)域畦田內(nèi)的灌水量之和。田間水利用系數(shù)是指根系層內(nèi)實際灌入的水量與進入毛渠水量的比值。

采用式（19）與（20）分別計算試驗區(qū)灌溉定額模擬值V及其與實測值間的相對誤差RE

（20）

式中V為灌溉定額模擬值，m3；V為灌溉定額實測值，m3；V為模擬的任意畦田的灌水量，m3；為試驗田總面積，m2；L為模擬的畦田長度，m；為模擬的田塊數(shù)量。

采用式（21）與（22）分別計算基于實測值與模擬值的單個田塊根系層的土壤平均含水率，

（22）

式中θ為基于實測值的根系層土壤平均質(zhì)量含水率，%；θ為基于模擬值的計劃濕潤層土壤平均質(zhì)量含水率，%；θ為模擬的0～100 cm范圍內(nèi)以1 cm為空間步長的各層土壤質(zhì)量含水率，%；1、2、3、4為實測的各層土壤質(zhì)量含水率，%。土壤根系層設(shè)為100 cm。

采用式（23）與（24）分別計算基于實測值或模擬值的田間水利用系數(shù)及兩者間的相對誤差RE：

（24）

式中V為灌溉后試驗區(qū)土壤根系層增加的總水量，m3；η為基于實測值計算的田間水利用系數(shù)；η為基于模擬值計算的田間水利用系數(shù)。

3.3 模型驗證

基于上述模型參數(shù)，畦長方向和土壤垂向的空間離散步長分別取值2.0和0.01 m、土壤與地表水流時間離散步長均取值1 s時，采用由畦灌地表水流-土壤水動力學耦合模型于Monte-Carlo抽樣建立的區(qū)域畦灌模擬方法，對3次試驗方案下區(qū)域畦灌過程進行模擬。對3次典型區(qū)域畦灌實例而言，基于確定性模擬方法FIM的灌水效率分別為68.78%、65.23%、67.17%，灌水均勻度分別為63.45%、62.14%、63.87%，基于所建隨機模擬方法SRM的灌水效率均值（標準差）分別為66.84%（8.45%）、67.23%（10.63%）、67.59%（9.21%）和63.23%（7.89%），灌水均勻度分別為61.58%（8.69%）、62.72%（8.36%）?？梢钥闯?，F(xiàn)IM與SRM之間的灌水效率和灌水均勻度之間存在一定的差異，而SRM通過標準差等指標可以有效描述區(qū)域畦田灌水質(zhì)量的內(nèi)在變異情況，F(xiàn)IM僅能采用固定指標值描述區(qū)域灌水質(zhì)量。

表4給出了不同模擬方法下灌溉定額模擬值與實測值之間的相對誤差與田間水利用系數(shù)模擬值與實測值間的相對誤差?？梢钥闯?，SRM的灌溉定額最小相對誤差為9.95%，最大相對誤差為12.23%，而已有模擬方法FIM的灌溉定額最小相對誤差為14.15%，最大相對誤差為16.78%。這表明構(gòu)建的模型具備優(yōu)良的模擬區(qū)域畦灌過程的能力。2類模擬方法的灌溉定額模擬值與實測值間的誤差，源于對土壤特性、畦田規(guī)格和地表糙率等灌溉參數(shù)區(qū)域空間變異性的捕捉精度有所差異。

表4 基于模擬與實測的灌溉定額及田間水利用系數(shù)相對誤差

注：FIM為確定性灌溉模型；SRM為隨機參數(shù)灌溉模型。

Note: FIM is deterministic parameter irrigation model, SRM is stochastic parameter irrigation model.

SRM模擬情況下3個實例的田間水利用系數(shù)相對誤差范圍為8.39%～10.21%，而FIM的田間水利用系數(shù)相對誤差范圍為13.87%～15.88%。這表明考慮畦灌參數(shù)空間變異性的畦灌隨機模擬方法能夠更好地評估區(qū)域畦田灌水質(zhì)量。2類模擬方法的模擬值與實測值間的誤差，可能是因為模擬時將田塊土壤特性空間變異性平均化處理，導致土壤水力特性參數(shù)值被過高或過低估計，或?qū)⑻锩嫖⒌匦螛藴什罴安诼氏禂?shù)概化處理，導致田面微地形和糙率系與實際情況有所偏差，從而影響對田間水利用系數(shù)的精確計算。

圖3給出了灌溉后畦田根系層的土壤含水率模擬（SRM）與實測均值在空間上的累積概率分布情況。可以看出，雖然實測田塊土壤平均含水率累積速度遲滯于模擬的土壤平均含水率累積速度，但模擬值與實測值的概率分布趨勢具有良好的一致性。土壤平均含水率模擬值與實測值累計概率曲線間存在的誤差，可能是因為模擬時對畦田微地形狀況、田面糙率系數(shù)、土壤水力特性空間變異、入滲維數(shù)及作物吸水等因素概化所導致的。

基于區(qū)域農(nóng)田耕作措施和種植作物相同，將量化田面微地形狀況的田面相對高程標準差視為一致，并將作物生長狀況等因素不以田塊尺度加以區(qū)分，視當次模擬試驗的田面糙率系數(shù)作為常數(shù)，但實際上各個田塊的田面相對高程標準差存在差異，田面糙率系數(shù)因為作物長勢客觀存在的差異及田面相對高程空間分布的差異而存在不同。土壤水力特性空間變異性在單個田塊及區(qū)域尺度上客觀存在，但模型模擬時，這種空間變異性采用離散后分配給各個虛擬田塊的方式加以表達，未能精確反映實際上單個田塊土壤水力特性空間變異性對地表水流推進過程及入滲過程的影響。另外，數(shù)值模擬時畦灌地表水流推進過程可視為從畦首到畦尾的一維地表水流運動過程，且單個田塊內(nèi)土壤入滲參數(shù)在點尺度上處處相等，但實際畦灌時地表水流運動間或存在著沿畦寬方向的曲線或漩渦運動過程，導致地表水流推進鋒不能同時到達相同畦田長度處。而作物根系通過影響土壤結(jié)構(gòu)和直接吸收水分2種形式改變土壤水分運移路徑，但模擬灌溉過程時忽略作物短期內(nèi)的吸水過程。因此，模擬與實測土壤平均含水率累積概率分布之間的誤差，是必然存在的。這種誤差只能在建立模擬方法時，通過合理簡化區(qū)域畦灌過程，建立精確合理的水流運動方程加以降低，在數(shù)值模擬時，通過輸入切合區(qū)域畦灌實際環(huán)境的參數(shù)值加以縮小，但由于灌溉過程及參數(shù)的概化，并不能完全消除。

通過對灌水效率、灌水均勻度、灌溉定額、灌后土壤平均含水率累積概率分布及田間水利用系數(shù)的模擬與實測值間的對比分析，發(fā)現(xiàn)所建的區(qū)域畦灌數(shù)值模擬方法具有精確的模擬能力，但由于對土壤特性、田面糙率系數(shù)、畦田微地形狀況、土壤入滲維數(shù)問題及作物吸水過程的概化，模型模擬結(jié)果與實測結(jié)果之間存在一定的誤差。這需要在下一步研究中加以改進和完善。

4 結(jié) 論

采用Monte-Carlo抽樣將具有空間變異性的區(qū)域灌溉參數(shù)離散表征為若干灌溉參數(shù)樣本，分別輸入給畦灌地表水流-土壤水動力學耦合模型，以模擬區(qū)域畦灌過程、評價區(qū)域畦田灌水質(zhì)量。該模擬方法有助于降低區(qū)域畦灌參數(shù)空間變異性的精確捕捉受到人為因素干擾的程度?；?次區(qū)域畦灌試驗的實測數(shù)據(jù)和2個對比的區(qū)域畦灌數(shù)值模擬方法，驗證了上述所建模擬方法的模擬效果。結(jié)果表明，區(qū)域畦田長度、入畦單寬流量、土壤砂粒、粉粒、黏粒、土壤容重、土壤有機質(zhì)含量均服從正態(tài)分布。確定性模擬方法與所建方法模擬的灌水效率和灌水均勻度具一定差異。與對比的模擬方法相比，所建模擬方法的模擬值與實測值間的灌溉定額和田間水利用系數(shù)相對誤差均較小，分別為9.95%～12.23%和8.39%～10.21%，而對比的確定性模擬方法下的模擬值與實測值間的灌溉定額和田間水利用系數(shù)相對誤差則分別為14.15%～16.78%和13.87%～15.88%，畦田土壤平均含水率累積趨勢表現(xiàn)出良好的一致性。所建模擬方法有效縮小了區(qū)域灌溉參數(shù)空間變異性未能很好描述所帶來的模擬誤差范圍較大等問題，為區(qū)域畦田灌溉過程及灌水質(zhì)量的精確評估、區(qū)域畦田灌溉優(yōu)化設(shè)計和管理提供了技術(shù)指導。但由于所建模擬方法對土壤特性、田面糙率系數(shù)、畦田微地形狀況、土壤入滲維數(shù)問題及作物吸水過程的概化，模擬與實測值之間存在一定的誤差，這需要在下一步研究工作中加以改進和完善。

[1] Furman A. Modeling coupled surface-subsurface flow processes: A review[J]. Soil Science Society of America, 2008, 7(2): 741－756.

[2] Abbasi F, Feyen J, Roth R L, et al. Water ?ow and solute transport in furrow-irrigated ?elds[J]. Irrigation Science, 2003, 22(2): 57－65.

[3] Abbasi F, ?im?nek J, van Genuchten M Th, et al. Overland water ?ow and solute transport: Model development and ?eld-data analysis[J]. Journal of Irrigation and Drainage Engineering, 2003, 129(2): 71－81.

[4] Zerihun D, Furman A, Warrick A W. Coupled surface-subsurface flow model for improved basin irrigation management[J]. Journal of Irrigation and Drainage Engineering, 2005, 131(2): 111－128.

[5] Bautista E, Zerihun D, Clemmens A J. Iterative coupling strategy for surface-subsurface flow calculations in surface irrigation[J]. Journal of Irrigation and Drainage Engineering, 2010, 136(10): 692－703.

[6] Banti M, Zissis Th, Anastasiadou-Partheniou E. Furrow irrigation advance simulation using a surface-subsurface Intera- ction model[J]. Journal of Irrigation and Drainage Engineering, 2011, 137(5): 304－314.

[7] 董勤各，許迪，章少輝，等. 一維畦灌地表水流-土壤水動力學耦合模型Ⅰ：建模[J]. 水利學報，2013，44(5)：570－577.

Dong Qin’ge, Xu Di, Zhang Shaohui, et al. Coupled surface-subsurface hydrodynamic model for border irrigationⅠ: Model development[J]. Journal of Hydraulic Engineering, 2013, 44(5): 570-577. (in Chinese with English abstract)

[8] 董勤各，許迪，章少輝，等. 一維畦灌地表水流-土壤水動力學耦合模型II：驗證[J]. 水利學報，2013，44(6)：726－733.

Dong Qin’ge, Xu Di, Zhang Shaohui, et al. Coupled surface-subsurface hydrodynamic model for border irrigation Ⅱ. Validation [J]. Journal of Hydraulic Engineering, 2013, 44(6): 726－733. (in Chinese with English abstract)

[9] 章少輝，許迪，李益農(nóng). 基于混合數(shù)值解法的二維畦灌地表水流運動模擬I：模型建立[J]. 水利學報，2011，42(2)：180－186.

Zhang Shaohui, Xu Di, Li Yinong. A 2D surface water flow simulation of basin irrigation based on a hybrid numerical method I. Establishment of model[J]. Journal of Hydraulic Engineering, 2011, 42(2): 180－186. (in Chinese with English abstract)

[10] Playan E, Slatni A, Castillo R, et al. A case study for irrigation modernization: II Scenario analysis[J]. Agricultural Water Management, 2000(42): 335－354.

[11] Paydar Z, Gallant J. A catchment framework for 1-D models introducing FLUSH and its application[J]. Hydrol Process, 2008, 22(13): 2094－2104.

[12] Awan U K, Tischbein B, Martius C. Combining hydrological modeling and GIS approaches to determine the spatial distribution of groundwater recharge in an arid irrigation scheme[J]. Irrigation Science, 2013, 31(4): 793－806.

[13] Skonard C J. A Field-scale Furrow Irrigation Model[D]. Nebraska: University of Nebraska, 2002.

[14] Minasny B, Mcbratney A. A conditioned Latin hypercube method for sampling in the presence of ancillary information[J].Computers & Geosciences, 2006, 32(9): 1378－1388.

[15] Rad M R P, Toomanian N, Khormali F, et al. Updating soil survey maps using random forest and conditioned Latin hypercube sampling in the loess derived soils of northern Iran[J]. Geoderma,2014, 232/234(12): 97－106.

[16] Shields M D, Zhang Jiaxin. The generalization of Latin hypercube sampling[J]. Reliability Engineering and System Safety, 2016, 148(12): 96－108.

[17] Goda T, Sato K. History matching with iterative Latin hypercube samplings and parameterization of reservoir heterogeneity[J]. Journal of Petroleum Science and Engineering, 2014, 114(2): 61－73.

[18] 章少輝，許迪，李益農(nóng). 基于混合數(shù)值解法的一維全水動力學模型[J]. 農(nóng)業(yè)工程學報，2009，25(9)：7－14.

Zhang Shaohui, Xu Di, Li Yinong. One-dimensional complete hydrodynamic model for border irrigation based on hybrid numerical method[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2009, 25(9): 7－14. (in Chinese with English abstract)

[19] García-Navarro P, Playán E, Zapata N. Solute transport modeling in overland flow applied to fertigation[J]. Journal of Irrigation and Drainage Engineering, 2000, 126(1): 33－40.

[20] Strelkoff T S, Clemmens A J, Bautista E. Field properties in surface irrigation management and design[J]. Hournal of Irrigation and Drainage Engineering, 2009, 135(5): 525－536.

[21] Mulder V L, de Bruin S, Schaepman M E. Representing major soil variability at regional scale by constrained Latin Hypercubeampling of remote sensing data[J]. Journal of Petroleum Science and Engineering, 2013, 21(1): 301－310.

[22] Bai Meijian, Xu Di, Li Yinong, et al. Stochastic modeling of basins microtopography: Analysis of spatialvariability and model testing[J]. Irrigation Science, 2010, 28(2): 157－172.

[23] 盛驟，謝式千. 概率論與數(shù)理統(tǒng)計[M]. 北京：高等教育出版社，1993.

[24] Schaap M G, Leij F J, van Genuchten M T. ROSETTA: A computer program for estimating soil hydraulic parameters with hierarchical pedotransfer functions[J]. Journal of Hydrology, 2001, 251(3/4): 163－176.

[25] Saxton K E, Rawls W J. Soil water characteristic estimates by texture and organic matter for hydrologic solutions[J]. Soil Science Society of America Journal, 2006, 70(5): 1569－1578.

[26] 白美健，李益農(nóng)，涂淑芳，等. 畦灌關(guān)口時間優(yōu)化改善灌水質(zhì)量分析[J]. 農(nóng)業(yè)工程學報，2016，32(2)：105－110.

Bai Meijian, Li Yinong, Tu Shufang, et al. Analysis on cutoff time optimization of border irrigation to improve irrigated water quality[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2016, 32(2): 105－110. (in Chinese with English abstract)

[27] 白美健，許迪，李益農(nóng). 微地形分布差異對畦灌過程及性能的影響模擬[J]. 農(nóng)業(yè)工程學報，2008，24(11)：1－6.

Bai Meijian, Xu Di, Li Yinong. Simulating effects of microtopography distribution difference on basin irrigation process and performances[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2008, 24(11): 1－6. (in Chinese with English abstract)

[28] 鮑士旦. 土壤農(nóng)化分析[M]. 北京：中國農(nóng)業(yè)出版社，1999.

[29] Lilliefors H W. On the Kolmogorov-Smirnov test for normality with mean and variance unknown[J]. Journal of the American Statistical Association, 1967, 62(318): 399－402.

[30] Bautista E, Clemmens A J, Strelkoff T S, et al. Modern analysis of surface irrigation systems with WinSRFR[J]. Agricultural Water Management, 2009a, 96(7): 1146－1154.

[31] Bautista E, Clemmens A J, Strelkoff T S, et al. Analysis of surface irrigation systems with WinSRFR: Example application[J]. Agricultural Water Management, 2009b, 96(7): 1162－1169.

Simulation and validation for regional border irrigation considering spatial variability of irrigation parameters

Dong Qin’ge1,2, Xu Di3※, Zhang Shaohui3, Bai Meijian3, Li Yinong3

(1.,712100,; 2.712100; 3.100038,)

The performance evaluation of regional scale border irrigation plays an important role in the improvement of surface irrigation management level, but the existing models present the shortcomings such as less accurately capturing the spatial variability of regional scale irrigation parameters. So the exiting models cannot be effectively applied to analyze the performance of regional scale border irrigation system. To solve this problem, the border surface-subsurface flow model was applied to describe the surface and subsurface water flow. The conservative complete hydrodynamic equation and the Richards equation were applied to describe the surface water flow and subsurface water flow in border irrigation, respectively. The finite-volume approach was applied to spatially discretize these governing equations to obtain good mass conservation ability. The Picard iteration approach was introduced to obtain the linearization of this nonlinear algebraic system. The mass conservation component of surface flow model and subsurface flow model were iteratively coupled at the same time step to obtain the convergence value of surface flow depth, and then the momentum conservation component of surface flow model was externally coupled based on the convergence value of both the surface flow depth and infiltration rate to update the surface flow velocity. Solutions were numerically computed using an improved hybrid numerical method for surface flow model and a proposed numerical solution method with high-order accuracy for subsurface flow model. Monte-Carlo sampling method was used to accurately capture the spatial variability of regional scale irrigation parameters and generate a large number of border irrigation parameters samples, which were input to border surface-subsurface model, respectively. Consequently, the border surface-subsurface water flow processes of regional scale could be accurately simulated.Three times border irrigation experiments at regional scale were performed to validate the proposed model in Mawan Irrigation District, located in Dongying City, Shandong Province. Soil samples were collected at 4 depths from the top, middle, and bottom of the border field to analyze soil bulk density, soil particle size distribution, and soil moisture. The soil hydraulic parameters were transformed from the abovementioned soil properties using the Rosetta model. The relative elevation values of the border bottom were observed using a water level gauge. The surface flow depth was measured using water depth measuring devices, which were placed at every observation point before the irrigation was initiated. The surface flow depth was used to estimate Manning’s roughness coefficient. And unit discharge, border length, and border width were observed in March 2012, November 2012, and March 2013. The validation results showed that the proposed model could well simulate regional scale border irrigation processes, and presented very good modeling accuracy. Specifically, the relative errors between the measured and simulated values by the proposed stochastic parameter irrigation model were 9.95%-12.23% and 8.39%-10.21% for irrigation quota and field water utilization coefficient, respectively. By contrast, the relative errors of irrigation quota and field water utilization coefficient based on the existing deterministic parameter irrigation models were 14.15%-16.78% and 13.87%-15.88%, respectively. Additionally, the cumulative probability distribution trends of average soil moisture after irrigation between the measured and simulated values present the satisfactory performance. Thus, the proposed model overcomes the shortcomings of existing models and provides a useful numerical analysis tool for the management and design of regional scale border irrigation system.

irrigation; Monte-Carlo method; soils; parameters; spatial variability; irrigation performance; numerical simulation

10.11975/j.issn.1002-6819.2017.09.001

S275.3; O212.2

A

1002-6819(2017)-09-0001-09

2016-07-22

2017-03-10

國家自然科學基金項目（51609237）；中科院“西部之光”項目（XAB2015B04）；國家863計劃項目（2013AA102904）；中央高校基本科研業(yè)務(wù)費專項資金項目（2452016104）

董勤各，男，河南南陽人，助理研究員，博士，主要從事畦溝灌節(jié)水技術(shù)、土壤水鹽運移模擬研究。楊凌 西北農(nóng)林科技大學水土保持研究所，712100。Email：qgdong2011@163.com

許迪，男，北京人，教授級高工，博士，主要從事灌溉技術(shù)理論與方法研究。北京 中國水利水電科學研究院，100038。 Email：xudi@iwhr.com

董勤各，許 迪，章少輝，白美健，李益農(nóng). 考慮灌溉參數(shù)空間變異的區(qū)域畦灌模擬與驗證[J]. 農(nóng)業(yè)工程學報，2017，33(9)：1－9. doi：10.11975/j.issn.1002-6819.2017.09.001 http://www.tcsae.org

Dong Qin’ge, Xu Di, Zhang Shaohui, Bai Meijian, Li Yinong. Simulation and validation for regional border irrigation considering spatial variability of irrigation parameters[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2017, 33(9): 1－9. (in Chinese with English abstract) doi：10.11975/j.issn.1002-6819.2017.09.001 http://www.tcsae.org