基于數值仿真的單目標多變量渠道控制參數尋優算法研究

張雨萌，管光華

（武漢大學 水資源與水電工程科學國家重點實驗室，武漢 430072）

0 引 言

【研究進展】自20 世紀30年代水力自動閘門的應用開創了渠道系統自動化控制的先河以來，國內外學者不斷研究創新，使得渠系自動化控制理論日趨成熟，相關的技術得到了一定范圍的應用。渠系自動化控制的核心是控制算法，即描述渠道輸入（水位或流量誤差）與輸出（閘門動作）的邏輯關系的方法。渠系控制算法先后出現了基于經典控制理論的渠道控制邏輯（如小人算法、科文算法以及后來的PID 類算法等[1]）；基于現代控制理論的渠道線性二次最優控制算法、最優下游水位控制器等；以及基于智能控制算法的魯棒算法等[2]。

近年來，渠系控制理論的各個方面又有了一些新的進展。姚雄等[3]對大型渠道控制蓄量的實現方法進行了研究，建立了包含倒虹吸、渡槽、取水口等建筑物在內的大型輸水渠道自動運行控制仿真模型，并對京石段進行了仿真驗證。Clemmens[4]提出了一種解決用水管理部門與用戶的供需不匹配問題的基于相鄰渠池的水位差異的控制器，其本質是將渠道系統看作一個大型水庫。Belaud 等[5]針對渠系響應時間展開研究，運用IDZ模型尋找基于用水需求改變的開環系統的恰當時間，并將其在一個渠系系統的不同的下游邊界條件下與SV 方程做比較，提出了灌溉渠系的響應時間概念。崔巍等[6]提出了一種具有較好實用性的閘前常水位運行控制算法，綜合前饋反饋算法提出動態調節渠道蓄量實現閘前常水位運行的控制算法，在南水北調的仿真測試中取得了良好效果。管光華等[7]對蓄量階躍補償、動力波原理、水量平衡模型3種時滯參數算法在ASCE的測試渠系和漳河灌區上進行仿真模擬，比較了它們的控制效果，從而推薦蓄量階躍補償算法為主要的前饋控制算法，并為簡化前饋規則提出了時滯參數顯式算法。

【擬解決的關鍵問題】上述控制理論及算法尚處于試驗和研究階段，在實際工程中，經典控制應用較為廣泛，特別是原理較為簡單、魯棒性較強的PID算法[2,8-12]。該算法在確定控制器參數時，主要有2種方法：一是傳統的經典設計方法，主要包括構建ID模型、IDZ模型等，并將其結合最優控制理論推算得到PID 參數；二是進行程序仿真試算，用圣維南方程結合尋優程序確定較優PID 參數。但由于方法一過于理想，沒有考慮真實渠道的局限性、非線性以及系統的通信延遲，故基于經典控制理論的設計方法、模型所得結果往往不夠理想，甚至出現控制失穩的情況。【切入點】目前，對于渠道控制系統PI 控制參數的尋優大多采用方法二中的網格法，本文采用仿真模型來率定和優化控制器參數，旨在將一些在其他領域運用較為成熟的算法（如梯度算法、PSO算法）改進運用于渠道控制系統中。【研究意義】以期能夠在渠系運行過程中較快地取得較好的PI 控制參數，使渠道達到控制目標，安全可靠及時地保障渠道系統高效運行。

1 PID 控制邏輯及其渠道控制性能指標

PID算法又叫作比例積分微分算法。此處的應用原理為：傳感器實時監控下游水位，將該水位與預定的常水位做比較，將二者之間的差記為誤差值e(t)。然后將該誤差值進行比例積分微分的線性組合計算，找到一個輸出的閘門開度u(t)作為控制器輸出變量[2]，計算式為：

式中：Kp為比例項因子；Ti為積分項因子；Td為微分項因子；e(t)為水位誤差值（m）；u(t)為閘門開度（m）。比例項是為了減少水位的誤差值；積分項是為了減少水位的穩態誤差；微分項是為了加快閘門的動作量，本文只考慮比例項和積分項[13]。

2 單渠池控制參數尋優算法

選取美國加利福尼亞州ASCE的小型測試渠道的單個渠池作為本文初步的尋優仿真對象，渠道參數見表1。流量變化工況設置為設計流量的15%～30%的階躍流量變化。

表1 單個渠池工況仿真渠道參數Table1 simulation parameters of single channel pool

2.1 網格法尋優及優化目標

網格法是PI 參數尋優最為經典的算法，其優點在于只要計算時間足夠長，取足夠密的網格進行尋優，就一定可以得到性能指標較好的一組控制參數。本文先以經典網格法進行尋優，將其作為一個比較對象與其他PI 參數尋優算法在運算效率、控制效果上做出比較，來評判其他算法是否有效或更優。

參考ASCE 提出的基于水位誤差變化幅度、渠系達到穩定所需時間、閘門動作量等關鍵性控制目標的性能指標，包括水位誤差平方積分（ISE）、最大絕對誤差（MAE）、穩態誤差（STE）、絕對值誤差積分（IAE）、絕對流量變化積分（IAQ）、絕對閘門開度積分（IAW等）[14]。其中，ISE是用于衡量水位波動平穩性能的指標；IAQ是描述系統在過渡過程中對流量的控制能力的指標；IAW是衡量閘門的啟閉幅度和頻率的指標。ISE越小，水位波動越小。初步探究時，本文認為若能控制水位波動在一個合理的范圍內，就能夠使系統達到一個較好的被控狀態。故采用較為傳統的ISE（水位誤差平方積分）作為尋優指標，計算式為：

式中：Δt為時間間隔；T為總的仿真時間；yt為當前時刻水位；ytarget為目標水位；ISE為水位誤差平方積分（此處為去量綱化，且多渠池工況下取平均值）。

綜合參考文獻[15]，本文中Kp、Ti的取值范圍為（0,100）。采取參數計算步長為2，對Kp、Ti進行（0,100）內的單渠池尋優，需計算50×50 共2 500個性能指標，共耗時約20 h（用于運算的電腦配置如下：Intel(R)Xeon(R)Gold5118 處理器，32GB 內存）。若再從中選取區間進行加密，需要約160 h 才能得到理想的控制參數，而本文2.4 節的計算結果顯示，運用本文提出的2種單渠池尋優方法，最多只需要8 h 就可以找到同樣控制效果的控制參數。可見網格法所付出的時間代價十分大，亟須改進。

2.2 復合梯度算法尋優

2.2.1 梯度法

梯度法是數學中較為經典的算法，廣泛應用于軍事、工業等方面的尋優。它利用負梯度方向是下降最快的方向的原理進行尋優，其迭代公式為：

式中：a0為步長；p(k)為目標值y對于變量的負導數；k為尋優次數[16]。

2.2.2 復合梯度法及其尋優效果

尋優分為2 部分，即隨機二分尋找梯度法初值以及梯度法尋優，稱為復合梯度法。由于不合適的控制器參數會使得系統失穩或性能指標顯著劣化，且梯度法尋優需要一個較好的初值。故認為梯度法控制參數的用于尋優的初始值對應的性能指標應優于沒有反饋控制（只有前饋控制）時的性能指標。

首先為了適應梯度法的初值敏感性，先使得PI參數（Kp、Ti）在一定范圍內（0～100）隨機生成n×n個值，并計算對應的性能指標值。若得到的最小的性能指標仍大于沒有反饋控制時系統對應的性能指標，則先數量級二分再區間二分，不斷縮小逼近最優范圍直到找到合適的梯度法初值。然后進行梯度法運算，以負梯度方向為最速下降方向為基礎，使得PI 參數對應的性能指標越來越好，直至達到預期性能指標（目前經仿真實驗驗證，預期性能指標取無反饋控制對應的性能指標的1/50～1/100，渠系控制效果良好），不斷逼近最優的PI 參數。

程序設計為當型循環，采取MPIEC（第i次的最小性能指標和i-1次最小性能指標的差）和MPIEL（第i-1次的最小性能指標和i-2次最小性能指標的差）作為循環條件，當二者都達到較小的程度（取預期性能指標精度的1/10），認為尋優進入了一個較小值的平緩區域，則停止循環，并將性能指標達到預期值作為結束循環的另一個條件。當遇到控制出錯的情況，例如閘門上游水位低于下游水位、上游閘門全關等情況。自出現控制出錯的時刻起，將性能指標擴充到整個尋優時間段并將對應性能指標擴大1 000倍（經試驗，每次遇到控制出錯的情況都會加長尋優時間，且此組性能指標必然不能用）。

使用梯度法原理進行PI 參數的迭代，計算式為：

式中：i為迭代次數；a0、b0為步長；IK、IT為ISE對Kp、Ti的偏導數的負值。運用式（4）和式（5）時，將其中a0、b0采用控制數量級的方法，使得加數的第2 項的數量級近似比加數第1 項的數量級少1個。經過多次驗算，該方法較其他方法更加有效穩定。

圖1 過閘流量過程線 Fig.1 Flow process line

圖2 閘門開度過程線Fig.2 Gate opening process line

圖3 下游水位誤差過程線Fig.3 Error line of downstream water level

表2 梯度法時間及循環次數統計Table2 Statistics of gradient method time and cycle times

2.3 復合PSO算法尋優

2.3.1 PSO算法 PSO算法，即粒子群算法，是基于鳥群覓食原理設計的[17]算法。先一次性撒下許多粒子，隨機生成N行D列的位置和速度數組。通過計算每個粒子的評價指標，可以篩選出此時的全局最優點。再進行M次迭代，對每個粒子的位置和速度進行一定規律的運算。該運算公式有許多不同的公式，綜合考慮尋優時間與尋優效果[18]。本文選取了較為基礎經典的一種計算方法，計算式為：

式中：v為粒子的速度；w為慣性權重；x為粒子所在的位置；i為粒子序號；c1、c2為比例因子；y為存儲x值的臨時存儲變量；pg為當前最優控制參數的存儲變量。

關鍵參數選取如下：N為初次撒下的點，本文認為初次撒下的點的效率較低，迭代計算的效率較高，故初選N=20；D為計算的維度，本例的研究變量僅有Kp、Ti共2個值，故D=2；w為慣性權重，認為粒子下一時刻的運動速度受上一時刻的運動速度所影響。此處取常規值0.6；M為迭代次數，是根據PSO算法的核心原理進行迭代計算，故取M=80；c1為粒子受自身最優值影響的比例因子；c2為粒子受群體最優值影響的比例因子。二者之和一般為4，目前取二者均為2[19]。

按此公式對每個粒子進行速度和位置的改變，并計算對應的性能指標來進行評判。不斷循環，直至找到滿足條件的性能指標，并篩選出對應的PI 參數，完成粒子群算法尋優[20]。

2.3.2 復合PSO算法及尋優效果

為提高尋優效率，參考2.2 節先使用隨機二分法生成N個尋優粒子（粒子含有位置：即Kp、Ti值；速度：下一步尋優的方向和步長），再利用PSO算法原理，使得N個粒子以不同的權重沿著上一時刻前進的方向、自身最優的方向、種群最優的方法進行前進循環，直到找到符合預期性能指標的PI 參數。

隨機二分法設計同梯度法，粒子群算法設計為：定義種群數（N=20）、粒子維數（D=2）、自身最優值學習因子（c1=2）、群體最優值學習因子（c2=2）、迭代次數（M=50）。先隨機生成N個二位數組x(i，n)（i=1～n），其物理意義為：隨機生成20 組Kp、Ti，并用x(i,1)儲存Kp的值；x(i,2)儲存Ti的值。再隨機生成N個二維數組速度v，以后期改變對應的第i組Kp、Ti值。

若該區域僅從空調風機盤管的角度考慮，可以采用以上公共走道和電梯廳采用的聯網型風機盤管的做法。本案例采用另外一種做法，借鑒客房客控RCU控制的方式，既可以實現遠程對風機盤管的控制，實現效果等同于聯網型風機盤管，又可以通過設置紅外探測器聯動公共衛生間內的照明開關，以實現節能效果。裙樓公共衛生間燈光和空調控制布置如圖6所示。

計算每組PI 參數對應的性能指標。并將當前20組Kp、Ti值參數存入到y矩陣里，并從中選擇對應性能指標最小的那組性能參數對應的PI 參數存入全局最佳性能指標對應PI 參數值pg里。對應的最佳性能指標值存入到c變量里。

再執行當型循環，循環條件為：迭代進行50次以內、性能指標達到預期值或者尋優進入平滑區域。多次根據粒子群算法核心公式更新鳥群的速度和位置。并不斷以更小的性能指標賦給變量c，不斷以對應的PI 參數賦給變量pg。若跳出當前循環后發現性能指標沒有達到預期值，則重新尋優。以此重復，直到鳥群的最佳性能指標小于預期值為止。最終的位置為pg，pg的第1 列為Kp，第2 列為Ti。c為最后的最佳性能指標[21-22]。

同樣，羅列出15 組計算結果如表3所示，表中TP表示粒子群法所用時間，NP表示粒子群法循環次數。水位誤差、閘門開度、過閘流量過程線顯示與復合梯度法得到的圖1—圖3控制效果相當：控制很快達到了穩定，水位誤差較小、控制效果也較為良好。

表3 PSO算法時間及循環次數統計Table3 Statistics of PSO algorithm time and cycle times

2.4 復合梯度法與復合PSO算法比較

本文以梯度算法與PSO算法尋優的各15 組算例進行相同尋優效果下尋優時間的對比，結果見表4。在控制效果相當的情況下，就其用時來說，2種算法較均勻網格法的效率均得到很大程度的提高，最多用時僅約為網格法的5%。但總體來說，復合梯度法所用平均時間僅僅為復合粒子群算法的35.76%，復合梯度法較復合粒子群算法更高效一些。二者對于初值的選取都較敏感，都需要先使用隨機二分法初步尋找主體算法初始值。二者的收斂性都較好，尋優一段時間均能找到良好的控制參數。

表4 2種算法尋優時間比較Table4 Comparison of optimization time between the two algorithms

3 多渠池工況下高效尋優

本文選取美國加利福尼亞州ASCE的小型測試渠道中連續3個渠池作為本文進一步的尋優仿真對象，渠道糙率為0.02，底坡為0.000 1，底寬為7 m，邊坡為1.5，渠道長度見表5。工況設置為30%～40%的階躍流量變化。將復合梯度法運用于多渠池工況，進行多變量尋優。并考慮到其復雜以及耗時性，先采用固定權重的多渠池多變量尋優，再采取隨機多渠池多變量尋優，由淺及深逐步探索。

3.1 固定權重的多渠池多變量尋優

將單渠池問題擴充到多渠池問題后，應該是6個變量的尋優，尋優難度會大大增加，效率也會迅速降低，故先采用簡化算法進行探索。即運算中第1個渠池的控制參數Kp、Ti仍然沿著負梯度方向進行尋優，而其他2個渠池的控制參數利用渠池之間的關系進行直接賦值（尋優指標為多渠池對應的性能指標）。并在其過程中，以單渠池梯度法為基礎，充分利用單渠池尋優所發現的結果和所利用的技巧，如解決梯度法初值敏感性的初值的隨機二分，選取步長前進所采用的控制數量級法，將肯定會出錯的性能指標根據時間占比擴充到整個時間段內等。

根據不同渠池的水面面積之間的關系決定不同渠池的控制參數的關系做法較為合理。根據相關文獻研究可知，渠池長度越長（水面面積越大），敏感性越小，對應的控制參數越大[23]。故本研究采用各個渠池的控制參數與其對應水面面積成正比的方式來選取不同渠池的Kp、Ti的值（Ti為反比關系）。

表5 多渠池各渠道長度 Table5 Length of channel

為了解決梯度法的初值敏感性問題，多渠池仍然采用先隨機二分，來取得大于無反饋控制時的性能指標所對應的尋優參數值，再將取得的初始值代入多渠池梯度法的核心算法中的方法，即隨機給定第1個渠池的Kp、Ti，再根據水面面積之比給定其他2個渠池的尋優參數。將這6個有一定關系的尋優參數先進行隨機二分法尋優，找到對應性能指標大于無反饋控制的一組值，作為梯度法的合理初始值代入多渠池的梯度算法，以期得到目標性能指標（此處的目標性能指標的數量級選取無反饋控制時性能指標的約1/100。若真正用于工程中，選取真實工況的預期性能指標精度即可）。

取得3個渠池的初始尋優參數后，沿著負梯度方向尋找第1個渠池的性能指標，而其他渠池的性能指標則采取每一步都根據水面面積之比給定。同單渠池算法，比較本次尋優參數所求得的性能指標和之前的性能指標，將更小值存在變量minJUG里。并求得MPISC、MPIEL不斷循環，直到MPISC和MPIEL都較小，說明尋優進入了一個較為平坦的區域停止循環。或者minJUG達到預期精度，則停止循環。最終輸出對應的性能指標和對應的尋優參數Kp1、Kp2、Kp3、Ti1、Ti2、Ti3。進入梯度法循環80次，若80次后仍未找到預期的性能指標，則進行重新尋優，重復以上步驟，不斷循環直到找到預期值。（根據實驗結果，此種情況循環2次左右即可找到預期值）完成多渠池梯度法運算。

得到的結果任取10 組取最大值、最小值、平均值列于表6。由表6可知，固定權重的多渠池多變量尋優法的平均尋優時間為4 h，較網格法大大提高了尋優效率。各項性能指標均較小、穩定時間較短、閘門運動幅度較小。從圖4—圖6所示過程線可得，水位誤差、流量過程變化均較為合理，控制效果良好。

3.2 隨機多渠池多變量尋優

為取得更好的控制效果，本文進行了6個變量的復雜多渠池尋優。仍然采用將尋優分成2個部分。第一部分為隨機二分粗選梯度法初值，第2 部分為梯度法尋優。與3.1 節的區別在于，第1 部分采取隨機生成6個毫無關系的變量，以及第2 部分梯度法尋優也是多渠池性能指標分別對6個變量求差分，6個變量都沿著自身的負梯度方向進行尋優。剩余算法處理同上一節大致相同，只是運算量由于求了3次差分而大大增加。

得到的結果任取10 組，其最大值、最小值、平均值列于表7。可見隨機多渠池多變量尋優性能指標較小、控制效果良好。只是尋優時間遠大于固定權重的多渠池多變量尋優。

圖4 過閘流量過程線 Fig.4 Flow process line

圖5 閘門開度過程線Fig.5 Process line of gate opening

圖6 水位誤差過程線Fig.6 Downstream water level error line

表6 固定權重梯度法尋優結果Table6 Optimization results of fixed weight gradient method

表7 多渠池多變量尋優結果Table7 Multivariate optimization results of multi-channel pools

表8 多渠池工況下2種算法結果比較Table8 Comparison of two algorithms under multi-channel pool conditions

3.3 2種算法比較

2種算法的各項性能指標結合表6—表8 進行比較，隨機多渠池多變量尋優的控制效果較固定權重的多渠池多變量尋優有較大的改善，隨機法較固定權重法的IAQ小90.48%、IAW小77.57%，Tstable（系統穩定時間）小17.05%，ISE小50.42%。根據2.1 節各性能指標所示物理意義可知，隨機法所得的控制器對流量的控制能力更強，閘門的起閉幅度和頻率更小，穩定時間更短，水位波動更小。但隨機法的所付出的時間代價也較大，較固定權重法多63.19%。用戶應綜合考慮利弊，根據自身需要來選取更為適用的方法。

4 討 論

目前，在水利工程領域研究PI 參數試錯法尋優的學者較少，其中，Stringam 等[12]提到的迭代法尋優就是實際應用中PI 參數尋優的經典方法。該法高度依賴經驗，需要渠道運行人員能夠依靠經驗給定PI參數較為合理的比例項因子初值。首先只進行比例控制，選取較為合理的比例項值；再固定比例項因子的值，以一定的規律改變積分項因子，直至渠系達到較好的控制效果。該法較為耗時，僅進行單渠池尋優也需要近8 h。對于缺乏良好管理經驗的管理人員的渠系，給定不恰當的PI 參數容易使渠系控制失調。且在較為緊急的渠系水量需求改變的工況下，也不能及時滿足用戶需求。而本文提出的幾種算法具有較好的普適性，對渠道運行管理人員經驗的要求不高，且控制的穩定性和安全性更有保障。

但本文的研究仍存在一些有待改進的方面：多渠池工況由于其存在不同渠池之間的耦合性，可能會造成誤差累計或放大。且隨著尋優渠道數目的增加，尋優變量也成倍增加，造成尋優時間大幅增加。故本文先采用固定權重的多渠池多變量尋優算法。之后，為尋找更好的控制指標，采用了隨機多渠池多變量尋優法，將所有變量均沿梯度法進行尋優，雖然得到的控制效果更好，但所付出的時間也多了約60%。在接下來的研究中，筆者將嘗試尋找控制參數與渠道特性及對應工況的更深層次的關系，以求能夠在盡量減少尋優時間的情況下，尋找控制效果更好的控制參數。

本文僅針對ISE進行尋優，只考慮了控制水位誤差較小這一優化目標。采用該性能指標的缺點是會造成閘門的過度調節，實際輸配水工程的運行調度需綜合考慮流量誤差、閘門運動、系統穩定時間等方面的要求。故可采用更加綜合的尋優指標，以獲得更為均衡的控制效果。

實際應用中，尋找仿真性能最優的控制器往往不是解決控制問題的良策，因糙率、閘門流量關系等水力要素的不確定性及時變性，“最優”的控制器在應用中并不一定有最佳的性能，而調度者更為關注的是控制性能在較大的范圍內能保持較好的水準，及控制器應具有較好的魯棒性。故可考慮進一步將控制性能的梯度引入到尋優目標中。

5 結 論

1）本文提出的復合梯度法和復合PSO算法都能較好地實現渠道控制系統PI 控制參數高效尋優。在控制效果相當的情況下，復合梯度法的尋優時間縮短為原來網格法的0.5%，PSO算法的尋優時間縮短為網格法的1.4%。

2）2種復合尋優算法都可以較好地解決尋優初值不穩定及尋優效率低下等問題。

3）多渠池工況下，隨機多渠池多變量尋優算法與固定權重的多渠池多變量尋優算法所取得的PI 控制參數都能較好地對渠道實現控制。隨機多渠池多變量尋優算法的控制效果明顯優于固定權重的多渠池多變量尋優算法。前者水位誤差平方積分（ISE）是后者的50%，穩定時間（Tstable）是后者的17%。但固定權重的多渠池多變量尋優算法的時間成本遠低于隨機多渠池多變量尋優，是隨機多渠池多變量尋優時間的37%。在進行控制器設計時，應視需求選擇合適的算法。