串列式雙風輪風電機組氣動特性建模及仿真

胡 陽, 王浩楠, 房 方, 許昱涵, 劉吉臻

(1.新能源電力系統(tǒng)國家重點實驗室,北京 102206;2.華北電力大學 控制與計算機工程學院,北京 102206)

隨著人類對能源需求的日益加劇,面對“碳達峰、碳中和”的戰(zhàn)略目標[1],建設(shè)新型清潔能源電力系統(tǒng)已是大勢所趨。風力發(fā)電(以下簡稱風電)是清潔能源的主要形式之一,風電機組具有建設(shè)周期短、機組壽命長、可靠性高和運行維護簡單的特點,因此受到廣泛關(guān)注[2]。風電機組中氣動系統(tǒng)通過降低空氣流速吸收空氣動能,將風能轉(zhuǎn)換為風輪旋轉(zhuǎn)的機械能[3];由于空氣動力學特性使風能捕獲過程出現(xiàn)復(fù)雜的本質(zhì)非線性[4],給風電過程的建模帶來了嚴峻的挑戰(zhàn),因此對風電機組氣動非線性系統(tǒng)的建模研究非常重要。

傳統(tǒng)風電機組氣動系統(tǒng)建模為葉素動量理論建模[5],該模型為白箱模型,結(jié)構(gòu)比較復(fù)雜。文獻[6]～文獻[7]研究了風電機組的氣動系統(tǒng)靜態(tài)非線性參數(shù)辨識問題,分別運用六參數(shù)模型和八參數(shù)模型進行研究,復(fù)雜度較高。氣動系統(tǒng)建模還有Cp-λ-β曲線建模,包括函數(shù)法和表格法,其中Cp為風能利用系數(shù),λ為葉尖速比,β為槳距角。該模型可以描述氣動系統(tǒng)的靜態(tài)特性,但無法表征2個靜態(tài)之間的超調(diào)等暫態(tài)現(xiàn)象。Zhao等[8]將風機非線性氣動模型識別并轉(zhuǎn)換為離散的基于聚類的線性分段仿射(PWA)模型,通過對多維工作域進行分類并優(yōu)化重組區(qū)域來估計線性子模型,該模型與某些選定工作點的線性化相比具有一定優(yōu)勢。潘晨陽[9]在建立氣動PWA模型表征氣動系統(tǒng)全局特性的同時,使用神經(jīng)網(wǎng)絡(luò)對誤差進行補償,結(jié)果表明使用該模型跟蹤氣動機械轉(zhuǎn)矩,可以改善模型精度。

新型串列式雙風輪風電機組在發(fā)電機的兩側(cè)裝設(shè)前、后2個風輪,相較于傳統(tǒng)單風輪機組,風能利用率有所提高。新型串列式雙風輪風電機組應(yīng)用了新型高效的對轉(zhuǎn)式異步發(fā)電機,結(jié)構(gòu)相對簡單,是一種具有廣闊應(yīng)用前景的新型風能轉(zhuǎn)換裝置。雙風輪風電機組氣動部分模型的搭建不同于單風輪風電機組,其前、后2個風輪之間存在無法忽略的氣動耦合關(guān)系[10],這使得雙風輪風電機組氣動部分具有高度復(fù)雜的非線性。因此,解決雙風輪風電機組氣動部分的非線性問題,對實現(xiàn)氣動部分的建模仿真具有重要意義。李新凱等[11]考慮前、后風輪各誘導因子和氣動耦合系數(shù),確定了雙風輪風電機組氣動系統(tǒng)性能與雙風輪風功率的變化趨勢,適用于雙風輪機組氣動性能設(shè)計。No等[12]利用雙風輪機組的一般運動方程,采用葉片單元和動量理論,以及流動相互作用原理構(gòu)造了氣動系統(tǒng)計算模型。史運濤等[13]推導了具有多輸入多輸出(MIMO)復(fù)雜形式的混雜系統(tǒng)的分段仿射自回歸(PWARX)模型的表達式并實現(xiàn)了模型辨識,該模型可應(yīng)用于各種具有復(fù)雜非線性特性的系統(tǒng),但目前鮮有相關(guān)研究提出雙風輪風電機組氣動部分PWARX建模表述。

筆者提出一種以雙風輪風電機組為研究對象的氣動MIMO系統(tǒng)全工況PWARX建模方案。該方案考慮輸入輸出數(shù)據(jù)的延遲影響,通過雙風輪風電機組機理分析結(jié)合系統(tǒng)辨識定階構(gòu)建有限差分回歸向量,針對雙風輪風電機組全工況運行數(shù)據(jù),以前、后風輪氣動轉(zhuǎn)矩為輸出,采用高維聚類直接對有限差分空間凸劃分并估計子模型邊界超平面,得到若干作用域。根據(jù)模型結(jié)構(gòu)辨識線性子模型參數(shù),集成實現(xiàn)雙風輪風電機組氣動MIMO系統(tǒng)全工況動態(tài)表征。最后,基于所搭建的雙風輪風電機組半物理仿真平臺進行仿真分析,以驗證該方案的有效性。

1 MIMO形式的PWARX建模結(jié)構(gòu)

以PWARX模型為基本結(jié)構(gòu),假設(shè)輸入量含h維,即u(k)=[u1,u2,…,uh]T∈Rh,輸出量含d維,即y(k)=[y1,y2,…,yd]T∈Rd,該形式單輸出的S個子模型結(jié)構(gòu)[14-15]可表示為:

(1)

式中:yg(k)為模型第g個輸出量,g=1,…,d;μg,i=[ag1,iT,ag2,iT,…,agna,iT,bg1,iT,bg2,iT,…,bgnb,iT,fg,i]T為模型參數(shù)向量,其中ag1,i,ag2,i,…,agna,i和bg1,i,bg2,i,…,bgnb,i分別為第i個子模型各階輸出和輸入對應(yīng)的系數(shù),維數(shù)為d維和h維,fg,i為偏移量;x(k)=[yT(k-1),…,yT(k-na),uT(k-1),…,uT(k-nb)]T=[x1(k),x2(k),…,xi(k),…,xn(k)]T為有限差分回歸向量,n=dna+hnb;na和nb分別為模型輸出量和輸入量的延遲階次;χi為第i個子模型作用域。

基于該模型結(jié)構(gòu)直接對多個輸出建模,得到的MIMO系統(tǒng)PWARX模型可描述為:

(2)

根據(jù)已知的數(shù)據(jù)集計算模型的參數(shù)向量和作用域,建模步驟如下:

(1) 考慮MIMO數(shù)據(jù)間延遲特性,定義含多個輸出的有限差分回歸向量x(k),找出每個數(shù)據(jù)點與附近擁有最小歐氏距離的c-1個點建立局部數(shù)據(jù)子集Ck。在每個數(shù)據(jù)子集中形成特征向量ξk,ξk由數(shù)據(jù)子集的參數(shù)向量λk結(jié)合數(shù)據(jù)子集中各點的均值θk構(gòu)成,即ξk=[λkT,θk]T。

基于文獻[15]中對λk的算式分析,推導出直接計算MIMO系統(tǒng)數(shù)據(jù)集參數(shù)向量的方法。該數(shù)據(jù)集參數(shù)向量的維數(shù)與輸出量個數(shù)有關(guān),通過最小二乘法計算,表達式如下:

(3)

式中:yCk,i為Ck中樣本第i個輸出量。

(2) 計算λk的經(jīng)驗協(xié)方差矩陣Vk[15],用各數(shù)據(jù)集回歸向量計算各類數(shù)據(jù)的離散度矩陣Lk:

(4)

將特征向量ξk看成遵循高斯分布的隨機向量,其協(xié)方差可估計為Rk=[Vk,0;0,Lk],使用K-means算法[16]對特征向量聚類,構(gòu)造S組有限差分數(shù)據(jù)集,記為D1,D2,…,DS。

(3) 把聚類得到的S組數(shù)據(jù)作為各子模型參數(shù)辨識時使用的數(shù)據(jù),以特征向量的置信度ωk[17]作為數(shù)據(jù)集的權(quán)重,采用加權(quán)最小二乘法辨識可得到各子模型的參數(shù)μi。子模型作用域的邊界可采用軟間隔支持向量機(SVM)求取超平面方程的系數(shù)[18]。算法優(yōu)化目標公式如下:

s.t.tk(φTxk+q)≥1-δk,δk≥0

(5)

式中:J為目標函數(shù);φ為不同工作域切平面的法向量;η為范圍可從0變化到1的懲罰系數(shù);δk為松弛變量,反映了數(shù)據(jù)不滿足硬間隔約束的水平;tk為可取不同數(shù)值的數(shù)據(jù)分類標志;q為偏移量;m為總數(shù)據(jù)量。

2 雙風輪機組氣動特性建模

將上述模型結(jié)構(gòu)應(yīng)用于雙風輪風電機組氣動特性,雙風輪氣動轉(zhuǎn)矩的輸出值與前后風輪輸入輸出數(shù)據(jù)具有復(fù)雜的非線性關(guān)系,可表示為:

Tr=fTr(V,ωr1,ωr2,β1,β2,Tr1,Tr2)

(6)

式中:Tr為模型氣動轉(zhuǎn)矩輸出值;Tr1、Tr2分別為前、后風輪氣動轉(zhuǎn)矩,N·m;V為來流風速,m/s;ωr1、ωr2分別為前、后風輪轉(zhuǎn)速,r/min;β1、β2分別為前、后風輪槳距角,(°);fTr為非線性函數(shù)。

考慮前、后風輪間的干涉特性,設(shè)置輸入向量含5維,u(k)=[u1,u2,u3,u4,u5]=[V,ωr1,ωr2,β1,β2]T∈R5;輸出變量含2維,y(k)=[y1,y2]=[Tr1,Tr2]T∈R2。根據(jù)上述建模方案,建立雙風輪氣動MIMO系統(tǒng)全工況PWARX模型,具體建模過程如圖1所示。

圖1 雙風輪機組氣動MIMO系統(tǒng)全工況PWARX建模方案

采用上述模型結(jié)構(gòu)定義有限差分回歸向量時,考慮各回歸向量中不同變量單位不同、數(shù)值相差很大等因素的影響,為了提高模型精確度,在建模前對回歸變量歸一化[19]處理表示如下:

(7)

式中:xnorm(i)為歸一化值;xmin(i)和xmax(i)分別為變量x(i)的最小值和最大值。

在進行高維聚類的過程中,當高維空間聚類取得的效果較好時,卻可能在三維空間展示時表現(xiàn)為不滿足預(yù)期的情況,此時可采用聚類性能評價指標量化,如果指標尚可,那么聚類效果就認為達到預(yù)期。引入戴維森堡丁指數(shù)(DBI),DBI指標RDBI考慮聚類數(shù)據(jù)簇內(nèi)緊密性和簇間分散性,越小表明聚類效果越好,定義如下:

(8)

式中:N為聚類簇個數(shù);Si、Sj分別為簇i、簇j數(shù)據(jù)各點到簇質(zhì)心的平均距離,用于度量樣本間的緊密程度;Mi,j為簇i與簇j質(zhì)心間的距離,用于度量不同簇樣本間的分散程度。

根據(jù)式(2)的結(jié)構(gòu)對各子模型參數(shù)進行辨識并劃分作用域后,集成子模型逼近雙風輪氣動全局非線性特性建立雙風輪氣動MIMO系統(tǒng)的PWARX模型。為評價模型性能引入模型性能評價指標,定義均方根誤差(ORMSE)、平均絕對誤差(OMAE)和平均絕對百分比誤差(OMAPE)。指標越小,表明模型性能越準確。

(9)

(10)

(11)

式中:ye(i)為預(yù)測值;y(i)為實際值。

3 仿真分析

3.1 半物理仿真平臺設(shè)計

為驗證所提方法的有效性,選取雙風輪風電機組半物理仿真平臺運行數(shù)據(jù)進行仿真驗證,平臺搭建一套不改變雙風輪運行原理及能量流動方向的軟硬件設(shè)備并設(shè)計相應(yīng)的控制系統(tǒng),仿真雙風輪機組發(fā)電過程,驗證其運行性能。平臺主要包括風機仿真模型、硬件可編程邏輯控制器(PLC)主控系統(tǒng)和上位機監(jiān)控系統(tǒng),各系統(tǒng)間使用OPC通信協(xié)議進行變量的通信交互,通信采樣間隔為0.1 s,設(shè)計架構(gòu)如圖2所示。

圖2 半物理仿真平臺架構(gòu)設(shè)計圖

半物理仿真平臺中的風機仿真模型設(shè)計首先基于雙風輪風電機組有限元分析建模仿真和尋優(yōu),得到雙風輪非線性氣動特性參數(shù),見表1。

表1 雙風輪風電機組氣動參數(shù)

利用前后風輪氣動參數(shù)、氣動干涉特性和葉素動量理論建立雙風輪葉素動量模型,基于該葉素動量模型設(shè)計流體力學特性仿真實驗,得到若干組雙風輪運行數(shù)據(jù),根據(jù)其中一組實驗數(shù)據(jù),得到風速范圍為3～15 m/s、精度為1 m/s時理想前后風輪轉(zhuǎn)速和槳距角,使機組最大程度利用風能獲得最大的風能利用系數(shù)(Cp,定義為Cp1+Cp2,其中Cp1、Cp2分別為前后風輪的風能利用系數(shù)),得到最佳運行工況,見表2。

表2 最佳運行工況表

根據(jù)表2,以及其他幾組不同物理量、維度及精細度的實驗數(shù)據(jù)分析雙風輪風電機組運行特點,制定雙風輪風電機組全工況覆蓋的運行區(qū)間。整體區(qū)間劃分結(jié)果主要分為最大功率跟蹤區(qū)、過渡區(qū)和恒功率區(qū),如圖3所示。

圖3 雙風輪機組全工況運行區(qū)間

根據(jù)雙風輪風電機組全工況運行區(qū)間劃分結(jié)果,在硬件PLC中設(shè)計相應(yīng)的控制算法,并在上位機監(jiān)控系統(tǒng)中實現(xiàn)對風機仿真模型和硬件PLC的變量交互通信,完成整個雙風輪風電機組半物理仿真平臺的搭建。

3.2 建模數(shù)據(jù)采集

選取在半物理仿真平臺運行的3～22 m/s全工況風速場景下生成的17 127組運行數(shù)據(jù)為建模數(shù)據(jù),包括來流風速、前后風輪轉(zhuǎn)速、前后風輪槳距角和前后風輪氣動轉(zhuǎn)矩,如圖4所示。

(a)

3.3 基于高維聚類的PWARX建模

根據(jù)雙風輪機理信息結(jié)合系統(tǒng)辨識確定延遲階次na=2,nb=1。對數(shù)據(jù)歸一化后構(gòu)造有限差分空間x(k)=[Tr1(k-1),Tr1(k-2),Tr2(k-1),Tr2(k-2),V(k-1),ωr1(k-1),ωr2(k-1),β1(k-1),β2(k-1)]T。聚類時設(shè)置局部數(shù)據(jù)子集個數(shù)c=50,根據(jù)機組運行區(qū)間特性選擇工作域個數(shù)為S=3。

對特征向量聚類后,選取每類數(shù)據(jù)集有限差分回歸向量第5、6、8組數(shù)據(jù)展示前風輪聚類結(jié)果,如圖5(a)所示,選取第5、7、9組數(shù)據(jù)展示對應(yīng)后風輪聚類結(jié)果,如圖5(b)所示。

(a) 前風輪

經(jīng)聚類后,由聚類量化指標計算可知RDBI為0.141 5,因此認為此次聚類結(jié)果達到預(yù)期目標。各子模型按照式(12)的模型結(jié)構(gòu),基于各工作域數(shù)據(jù)子集對模型參數(shù)進行辨識,辨識結(jié)果見表3。針對相臨數(shù)據(jù)集采用軟間隔SVM得到邊界超平面方程ωTx+b=0,系數(shù)ωT和b的結(jié)果見表4。χ1和χ2分別為第1、2個子模型和第2、3個子模型的超平面系數(shù);根據(jù)上述子模型參數(shù)和邊界超平面對應(yīng)得到包含3個工作域的雙風輪風電機組氣動MIMO系統(tǒng)全工況PWARX模型。

表3 模型辨識各項系數(shù)

表4 子模型邊界超平面

(12)

3.4 模型驗證及性能評價

為驗證所提PWARX模型精確度,采用所搭建的雙風輪風電機組半物理仿真平臺進行分析驗證,根據(jù)平臺生成不同工作域場景的運行數(shù)據(jù)來驗證該模型的準確性和辨識精度,并對模型的誤差進行評價。場景一、二分別為低、高風速2種單工作域運行場景,用于驗證該模型單工作域子模型的預(yù)測效果;場景三為多工作域切換運行場景。各工作域場景生成的風速序列和PWARX模型在各場景下預(yù)測輸出值和實際輸出值的對比結(jié)果如圖6～圖14所示。

圖6 場景一風速序列

(a) 前風輪

(a) 前風輪

圖9 場景二風速序列

(a) 前風輪

(a) 前風輪

圖12 場景三風速序列

(a) 前風輪

(a) 前風輪

單工作域場景模型性能指標見表5。由表5可知,模型應(yīng)用在2種單工作域場景的前后風輪各性能指標相近,由于僅選取各單工作域1 000組數(shù)據(jù),模型精度較高。

表5 單工作域模型性能指標

多工作域切換運行場景模型跟蹤輸出的精度較高,并且能克服工作域切換的擾動,以此類推達到多工作域子模型集成實現(xiàn)全工況動態(tài)逼近。多工作域模型評價指標見表6。由表6可知,模型在多工作域切換狀態(tài)下也能保持良好的逼近性能,對實際輸出有較好的跟蹤性能。

表6 多工作域模型性能指標

綜上所述,結(jié)合不同場景下的仿真驗證情況,采用的雙風輪風電機組氣動MIMO系統(tǒng)全工況PWARX建模方案在不同場景下對輸出實際值的擬合程度都沒有出現(xiàn)較大的偏差,在不同場景下模型的輸出值具有較好的逼近能力,都能很好地跟蹤實際輸出值的變化趨勢,辨識參數(shù)準確度較高,對雙風輪機組氣動系統(tǒng)具有較好的擬合效果。

4 結(jié) 論

(1) 基于雙風機氣動子系統(tǒng)MIMO特性,提出基于高維空間聚類、超平面估計、有限差分作用域等概念建立雙風輪風電機組氣動MIMO系統(tǒng)全工況PWARX模型來估計實際輸出氣動轉(zhuǎn)矩,表征雙風輪氣動系統(tǒng)的全局動態(tài)特性。

(2) 基于平臺運行數(shù)據(jù)定義有限差分空間,劃分作用域,適用于氣動系統(tǒng)非線性問題。根據(jù)雙風輪風電機組運行區(qū)間特性定義聚類個數(shù)相較于人為劃分得到的聚類效果更好,工作域劃分結(jié)果比較合理。

(3) 通過半物理仿真平臺的雙風輪風電機組模型設(shè)計實驗和仿真驗證,針對不同場景驗證模型參數(shù)辨識結(jié)果的準確性。仿真結(jié)果表明,所提PWARX模型具有良好的區(qū)域動態(tài)表征能力,且其輸出氣動轉(zhuǎn)矩能夠較好地跟蹤實際輸出值,從而表明該方法的有效性。