陳 上，陳金保，王安林，丘濤基，陳齊燈，肖志懷

（1.武漢大學(xué)動(dòng)力與機(jī)械學(xué)院，武漢430072；2.福建華電金湖電力有限公司，福建三明353399）

0 引 言

隨著大規(guī)模風(fēng)能、太陽能等新能源電能的接入，輸電電網(wǎng)架構(gòu)和區(qū)域電網(wǎng)間互聯(lián)方式日趨復(fù)雜，在電力系統(tǒng)的規(guī)劃、調(diào)峰、安全備用、頻率電壓穩(wěn)定等方面帶來更多新的問題。作為電網(wǎng)中起主要調(diào)節(jié)作用的水電站，對其調(diào)節(jié)性能提出了更高要求。由于水輪機(jī)內(nèi)部流體運(yùn)動(dòng)復(fù)雜，各運(yùn)行變量間的動(dòng)態(tài)關(guān)系復(fù)雜、耦合程度高，難以用簡單的數(shù)學(xué)模型描述，目前只能通過真機(jī)試驗(yàn)、模型試驗(yàn)和近似解析方法來建立水輪機(jī)模型。在工程上通常采用穩(wěn)態(tài)特性近似代替動(dòng)態(tài)特性建立水輪機(jī)模型［1］，進(jìn)而對水輪機(jī)的穩(wěn)態(tài)特性、動(dòng)態(tài)特性、控制策略等進(jìn)行研究。水輪機(jī)建模主要包括外特性法［2-6］和內(nèi)特性法［7-9］，內(nèi)特性建模由于過度依賴水輪機(jī)幾何參數(shù)和結(jié)構(gòu)參數(shù)，且實(shí)際建模時(shí)很難獲得具體準(zhǔn)確的參數(shù)，故很少采用。外特性建模主要基于模型運(yùn)轉(zhuǎn)綜合特性曲線（或原型運(yùn)轉(zhuǎn)綜合特性曲線），位于模型綜合特性曲線內(nèi)的工況通?？奢^好表達(dá)，但小開度區(qū)間建模還需要進(jìn)一步深入研究；此外，也有學(xué)者提出其他水輪機(jī)建模方法，如適用于電力系統(tǒng)仿真的遞推水輪機(jī)線性模型［10］；反映系統(tǒng)能量特性的非線性水輪機(jī)哈密爾頓模型［11］。為了分析各類模型之間的建模精度與差異，為水輪機(jī)建模提供理論指導(dǎo)，本文結(jié)合目前常用的分段線性化模型、理想水輪機(jī)模型、神經(jīng)網(wǎng)絡(luò)模型和插值法模型，借助matlab 仿真平臺(tái)，針對實(shí)際機(jī)組開展仿真建模研究，對比分析不同建模方法在穩(wěn)態(tài)特性、動(dòng)態(tài)特性、建模復(fù)雜程度、模型仿真時(shí)間等方面效果與差異。

1 水輪機(jī)建模原理

1.1 分段線性化模型

水輪機(jī)動(dòng)態(tài)特性通常使用水輪機(jī)穩(wěn)態(tài)工況下的力矩與流量特性表示，其非線性穩(wěn)態(tài)表達(dá)式為式（1）和式（2）：

式中：Mt為水輪機(jī)力矩；a為導(dǎo)葉開度；H為水輪機(jī)水頭；z為槳葉角度；n為機(jī)組轉(zhuǎn)速；Q為水輪機(jī)流量。

將這兩個(gè)式子展開為泰勒級數(shù)并略去二次及以上高次項(xiàng)得式（3）與（4）。

為方便描述，將其偏導(dǎo)數(shù)替換成傳遞系數(shù)，如式（5）：

該方程描述穩(wěn)態(tài)工況附近微小偏差時(shí)力矩和流量的變化。即該方程的應(yīng)用條件是已知該工況的傳遞系數(shù)，且波動(dòng)幅度在工程允許范圍之內(nèi)。當(dāng)工況波動(dòng)過大時(shí)，該方程計(jì)算數(shù)據(jù)會(huì)與實(shí)際結(jié)果產(chǎn)生較大的誤差，故當(dāng)波動(dòng)較大時(shí)，此方法應(yīng)用受到限制。

分段線性化建模中，關(guān)鍵是工況的劃分和傳遞系數(shù)的計(jì)算。工況劃分通常以水頭和導(dǎo)葉開度等比例選取，工況劃分越細(xì)，模型精度越高。

傳遞系數(shù)計(jì)算通常采用兩點(diǎn)法或者曲線擬合法，適合計(jì)算有限工況點(diǎn)的傳遞系數(shù)，但計(jì)算繁瑣，難以清晰全面認(rèn)識(shí)傳遞系數(shù)隨工況變化規(guī)律，故而出現(xiàn)通過模型綜合特性曲線結(jié)合BP神經(jīng)網(wǎng)絡(luò)對數(shù)據(jù)處理后直接進(jìn)行求偏導(dǎo)數(shù)的處理方法，此方法計(jì)算既保證了計(jì)算精度又減少了計(jì)算量［12］。

計(jì)算傳遞系數(shù)，得出如圖1 所示的水輪機(jī)分段線性化模型方框圖［13］，圖1中Tw為水流慣性時(shí)間常數(shù)。

1.2 理想水輪機(jī)模型

主要分析調(diào)節(jié)系統(tǒng)其他環(huán)節(jié)時(shí)，一般采用理想水輪機(jī)模型表示水輪機(jī)特性。理想水輪機(jī)模型是對1.1 節(jié)分段線性化模型進(jìn)行簡化，應(yīng)用時(shí)不考慮槳葉角度，通過槳葉角度與導(dǎo)葉開度協(xié)聯(lián)，將轉(zhuǎn)漿式水輪機(jī)視為混流式水輪機(jī)，考慮到轉(zhuǎn)速對力矩和流量影響較小，理想水輪機(jī)模型中傳遞系數(shù)ex和eqx取為0，在額定工況附近，取ey=1.0，eh=1.5，eqy=1.0，eqh=0.5，將傳遞系數(shù)合并得到傳遞函數(shù)如圖2所示。

圖2 理想水輪機(jī)模型Fig.2 Ideal hydraulic turbine model

1.3 神經(jīng)網(wǎng)絡(luò)模型

隨著人工智能發(fā)展，基于神經(jīng)網(wǎng)絡(luò)技術(shù)的非線性建模技術(shù)得到廣泛應(yīng)用。針對BP 神經(jīng)網(wǎng)絡(luò)，只要有足夠多隱含層和隱含層節(jié)點(diǎn)，就可以逼近任意非線性映射關(guān)系，許多學(xué)者將神經(jīng)網(wǎng)絡(luò)引入水輪機(jī)建模之中［14-17］。

實(shí)時(shí)仿真中，利用神經(jīng)網(wǎng)絡(luò)實(shí)現(xiàn)水輪機(jī)非線性特性計(jì)算具有計(jì)算簡單、精度可控且導(dǎo)數(shù)連續(xù)等優(yōu)點(diǎn)。可供選擇的神經(jīng)網(wǎng)絡(luò)很多，如前饋神經(jīng)網(wǎng)絡(luò)、循環(huán)網(wǎng)絡(luò)和對稱連接網(wǎng)絡(luò)等，水輪機(jī)建模中使用前饋神經(jīng)網(wǎng)絡(luò)中BP 神經(jīng)網(wǎng)絡(luò)較多，其網(wǎng)絡(luò)結(jié)構(gòu)如圖3所示。

圖3 BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)圖Fig.3 BP neural network structure diagram

該網(wǎng)絡(luò)包含三輸入、單輸出和一個(gè)隱層，實(shí)現(xiàn)軸流轉(zhuǎn)漿式機(jī)組流量或力矩模擬，其隱層神經(jīng)元個(gè)數(shù)可調(diào)整，輸出層為線性，隱層神經(jīng)元傳遞函數(shù)為log-signoid 函數(shù)。神經(jīng)網(wǎng)絡(luò)算法表達(dá)式如式（6）：

式中：NNout為神經(jīng)網(wǎng)絡(luò)輸出；N11為單位轉(zhuǎn)速，w1i、w2i、w3i為隱層神經(jīng)元權(quán)值；bi為隱層神經(jīng)元閾值；w0i為輸出層權(quán)值；b0為輸出層閾值；W1、W2為權(quán)值向量。

通過神經(jīng)網(wǎng)絡(luò)學(xué)習(xí)可得到網(wǎng)絡(luò)相應(yīng)權(quán)值與閾值，從而建立基于神經(jīng)網(wǎng)絡(luò)的水輪機(jī)模型，實(shí)際應(yīng)用中，只需輸入單位轉(zhuǎn)速、導(dǎo)葉開度和槳葉角度，利用已建立的水輪機(jī)神經(jīng)網(wǎng)絡(luò)模型，就可計(jì)算水輪機(jī)輸出流量和力矩。水輪機(jī)神經(jīng)網(wǎng)絡(luò)模型如圖4所示。

圖4 神經(jīng)網(wǎng)絡(luò)模型方框圖Fig.4 Neural network model block diagram

神經(jīng)網(wǎng)絡(luò)樣本數(shù)據(jù)獲取方法如下：

神經(jīng)網(wǎng)絡(luò)建模的依據(jù)為機(jī)組模型綜合特性曲線，引入單位流量、單位轉(zhuǎn)速和單位力矩如式：

式中：M11為單位力矩；D為水輪機(jī)直徑；H為水輪機(jī)水頭；M為水輪機(jī)力矩；Q11為單位流量；Q為水輪機(jī)流量；n11為單位轉(zhuǎn)速；n為水輪機(jī)轉(zhuǎn)速。

單位轉(zhuǎn)速n11、單位流量Q11、導(dǎo)葉開度α、槳葉角度β以及效率η是神經(jīng)網(wǎng)絡(luò)訓(xùn)練的輸入數(shù)據(jù)。首先以n11和α為輸入，單位流量Q11為輸出得到第一組數(shù)據(jù)，訓(xùn)練出的神經(jīng)網(wǎng)絡(luò)即Qnet。其次以n11和Q11為輸入，效率η為輸出得到第二組數(shù)據(jù)。以一個(gè)擬合函數(shù)或者BP神經(jīng)網(wǎng)絡(luò)，使用第二組數(shù)據(jù)建立輸入n11、Q11與輸出η的關(guān)系式（10）。坐標(biāo)變換后進(jìn)行神經(jīng)網(wǎng)絡(luò)訓(xùn)練。利用關(guān)系式（10）與關(guān)系式（11），計(jì)算得出輸入導(dǎo)葉開度α與單位轉(zhuǎn)速n11對應(yīng)輸出單位力矩M11，即為神經(jīng)網(wǎng)絡(luò)Mnet。

式中：η為水輪機(jī)效率；α為導(dǎo)葉開度；γ=ρ g，γ為容重。

1.4 插值法水輪機(jī)模型

插值法水輪機(jī)模型與神經(jīng)網(wǎng)絡(luò)類似，也可描述水輪機(jī)非線性特性。神經(jīng)網(wǎng)絡(luò)模型是通過數(shù)據(jù)處理訓(xùn)練得到權(quán)值與閾值，在仿真時(shí)通過歸一化后的單位轉(zhuǎn)速n11與導(dǎo)葉開度α進(jìn)行矩陣計(jì)算，再反歸一化得出輸出單位力矩與單位流量。插值法模型通過模型綜合特性曲線數(shù)據(jù)計(jì)算得到插值表，在仿真時(shí)輸入單位轉(zhuǎn)速n11與導(dǎo)葉開度α，通過二元三次插值計(jì)算輸出力矩與流量。神經(jīng)網(wǎng)絡(luò)模型需要前期樣本數(shù)據(jù)多，但實(shí)時(shí)仿真輸出數(shù)據(jù)誤差相對小，插值法使用的插值表劃分以及不同插值方法都會(huì)對最終結(jié)果產(chǎn)生影響。

插值法水輪機(jī)模型同樣基于模型綜合特性曲線建模，按照水頭和開度劃分工況，計(jì)算出對應(yīng)工況點(diǎn)流量和力矩。仿真時(shí)輸入單位轉(zhuǎn)速與開度，采用二元三次插值法進(jìn)行插值計(jì)算出全工況的流量與力矩。

二元三次插值公式為式（12）。

其中k≠i，l≠j，xp＜x＜xp+1，yq＜y＜yq+1。

其原理如圖5所示。具體步驟為：

圖5 二元三次插值示意圖Fig.5 Binary cubic interpolation schematic diagram

設(shè)插值點(diǎn)為S，其坐標(biāo)為S（xs，ys，zs），

（1）在插值點(diǎn)S周圍找出9 個(gè)點(diǎn)A、B、C、D、E、F、G、H、I，且滿足以下關(guān)系式：

其中：xp＜xs＜xp+1，yq＜ys＜yq+1。

（2）分別對（A、D、G），（B、E、F），（C、F、I）進(jìn)行一元三點(diǎn)拉格朗日插值，得到點(diǎn)（U、V、W）。

（3）對點(diǎn)（U、V、W）再進(jìn)行一次一元三點(diǎn)拉格朗日插值，得點(diǎn)S，S點(diǎn)的值zs，s即為該插值點(diǎn)的函數(shù)值。

仿真時(shí)插值法模型與神經(jīng)網(wǎng)絡(luò)模型相似，都是輸入單位轉(zhuǎn)速與導(dǎo)葉開度，神經(jīng)網(wǎng)絡(luò)模型通過s 函數(shù)里編寫的歸一化及矩陣計(jì)算代碼計(jì)算得出輸出力矩與流量，插值法模型通過simulink 中2-D lookup table 模塊實(shí)現(xiàn)插值，輸出力矩與流量。模型方框圖如圖6所示。

圖6 插值法模型方框圖Fig.6 Block diagram of interpolation model

2 仿真及對比實(shí)驗(yàn)

2.1 模型參數(shù)

某機(jī)組為混流式機(jī)組，該機(jī)組各參數(shù)如表1，神經(jīng)網(wǎng)絡(luò)權(quán)值與閾值如表2，神經(jīng)網(wǎng)絡(luò)輸入輸出關(guān)系曲面如圖7 所示，插值表輸入輸出關(guān)系如圖8所示。圖中X軸為單位轉(zhuǎn)速，Y軸為導(dǎo)葉開度，Z軸為輸出流量或力矩。

圖7 流量與力矩神經(jīng)網(wǎng)絡(luò)Fig.7 Flow and torque neural network

圖8 流量與力矩插值表Fig.8 Flow and torque interpolation table

表1 額定工況參數(shù)Tab.1 Rated condition parameter

表2 神經(jīng)網(wǎng)絡(luò)權(quán)值與閾值Tab.2 Neural network weights and thresholds

調(diào)速器為并聯(lián)型調(diào)速器，PID 參數(shù)為KP=0.20，KI=0.05，KD=0.05。分別以功率上調(diào)1%和10%為階躍輸入得到仿真力矩輸出如圖9（a）所示，流量輸出如圖9（b）所示，水頭輸出如圖9（c）所示。

2.2 仿真結(jié)果分析

圖9（a）為模型力矩輸出，如圖所示其中4 種模型穩(wěn)態(tài)力矩相近，只有動(dòng)態(tài)過程存在一定區(qū)別，其中分段線性化模型動(dòng)態(tài)過程最快，理想水輪機(jī)模型動(dòng)態(tài)過程最慢。圖9（b）為模型流量輸出，如圖所示4種模型流量波動(dòng)基本不同，其中理想水輪機(jī)模型的波動(dòng)量最小，而神經(jīng)網(wǎng)絡(luò)模型和插值法模型的波動(dòng)量最大，且穩(wěn)態(tài)值相近。圖9（c）為模型水頭輸出，如圖9所示除插值法外水頭波動(dòng)穩(wěn)態(tài)后基本相同，插值法在波動(dòng)開始時(shí)水頭產(chǎn)生較大瞬時(shí)變化然后趨于平緩。為更便于觀察各參數(shù)變化，引入公式（13），對力矩、流量和水頭變化進(jìn)行對時(shí)間求導(dǎo)后，除以功率變化，處理得圖10和11。

圖9 仿真輸出Fig.9 The simulation output

式中：t為時(shí)間；ΔP為機(jī)組功率變化；φ1為力矩導(dǎo)數(shù)對功率變化相對值；φ2為流量導(dǎo)數(shù)對功率變化相對值；φ3為水頭導(dǎo)數(shù)對功率變化相對值。

由圖10 可知在功率波動(dòng)初期，插值法水輪機(jī)模型的φ1、φ2與φ3經(jīng)過較大波動(dòng)然后趨于平緩，其他3 個(gè)模型的φ1與φ2參數(shù)相近，只有φ3存在較大差別。由不同功率波動(dòng)對比可知：功率波動(dòng)對于φ波形變化影響較小。由圖11可知：分段線性化模型仿真時(shí)間約為0.22 s，理想水輪機(jī)模型為0.21 s，插值法模型接近0.53 s，約為前兩種方法二點(diǎn)五倍，神經(jīng)網(wǎng)絡(luò)模型仿真時(shí)間更長，為18 s左右。

圖10 不同功率波動(dòng)時(shí)φ值變化Fig.10 The φ value changes with different powerfluctuations

圖11 實(shí)際仿真時(shí)間Fig.11 Actual simulation time

3 結(jié) 語

本文建立了4種水輪機(jī)模型，并進(jìn)行了仿真對比分析，結(jié)果表明：4種不同模型都能有效仿真力矩穩(wěn)態(tài)特性，動(dòng)態(tài)過程相差不大，而流量仿真時(shí)4 種模型輸出相差較大。進(jìn)行水頭仿真時(shí)相比于其他3種方法，插值法模型水頭波動(dòng)會(huì)偏大，水頭波動(dòng)反映到蝸殼水壓后會(huì)造成仿真蝸殼水壓大于實(shí)際水壓，進(jìn)而使調(diào)節(jié)保證計(jì)算產(chǎn)生一定誤差。仿真穩(wěn)態(tài)力矩時(shí)四種模型精度相近，在動(dòng)態(tài)過程存在一定差別，在波動(dòng)瞬間若超出插值表邊界處會(huì)產(chǎn)生水頭瞬時(shí)較大波動(dòng)，產(chǎn)生一定誤差。

分段線性化模型在需要仿真時(shí)間較短和建模更簡單時(shí)是較好的選擇，當(dāng)更多考慮動(dòng)態(tài)過程精確度時(shí)，選擇神經(jīng)網(wǎng)絡(luò)模型會(huì)更好。 □