柏文超，劉穎明，王曉東，高 興，張書源

（沈陽工業(yè)大學(xué) 電氣工程學(xué)院，遼寧 沈陽 110870）

0 引言

在“雙碳”政策背景下，風(fēng)電行業(yè)的發(fā)展模式正逐步向精細(xì)化、集約化發(fā)展模式轉(zhuǎn)變，行業(yè)研究關(guān)注點已經(jīng)集中到增功提效方面[1]。

當(dāng)風(fēng)電機(jī)組運行在額定風(fēng)速以下時，通常認(rèn)為葉片方位位于0°是最優(yōu)槳距角位置；此時，風(fēng)能利用系數(shù)達(dá)到最大，輸出功率最大[2]。然而，受葉片安裝誤差、機(jī)組外界運行環(huán)境變化等因素影響，機(jī)組運行時會發(fā)生最優(yōu)槳距角與理論槳距角不一致的情況，從而造成風(fēng)能利用率的下降；所以，需要對最優(yōu)槳距角的位置（可正可負(fù)）進(jìn)行研究[3]。

辨識最優(yōu)槳距角，對提升風(fēng)電機(jī)組發(fā)電量具有一定的價值和意義。

為快速找到風(fēng)電機(jī)組最優(yōu)槳距角，文獻(xiàn)[4]提出了一種自尋優(yōu)算法。由于用該方法構(gòu)造的槳距角評價函數(shù)不具有穩(wěn)定性，所以當(dāng)風(fēng)速波動劇烈時，得到的結(jié)果誤差較大。

文獻(xiàn)[5]考慮空氣密度因素，提出了一種通過調(diào)整風(fēng)電機(jī)組安裝角的方法來增大機(jī)組輸出功率。該方法的應(yīng)用僅限于定槳距風(fēng)電機(jī)組，且未考慮額定風(fēng)速以下機(jī)組發(fā)電量。

文獻(xiàn)[6]針對葉片安裝角存在人為誤差的問題，通過設(shè)置不同的安裝角來對其進(jìn)行模擬，并以此研究了安裝角對風(fēng)電機(jī)組性能的影響。但是，文中只分析了安裝角誤差對機(jī)組性能的影響，并未提出針對性解決方法。

為達(dá)到提升機(jī)組發(fā)電量的目的，本文提出一種基于運行數(shù)據(jù)驅(qū)動的風(fēng)電機(jī)組最優(yōu)槳距角辨識方法：首先，用擬合方法得到機(jī)組“風(fēng)速-功率-槳距角”連續(xù)三維特性曲面，然后求解三維特性曲面變量之間的函數(shù)關(guān)系，進(jìn)而辨識出不同平均風(fēng)速下風(fēng)電機(jī)組最優(yōu)槳距角。

1 數(shù)據(jù)驅(qū)動辨識框架

風(fēng)電機(jī)組槳距角對輸出功率的影響如圖1所示。從圖1可以看出：在同一風(fēng)速下，不同槳距角對應(yīng)的輸出功率不同。所以對于不同工況，可以通過改變機(jī)組槳距角的大小提升發(fā)電量。

圖1 槳距角與功率的關(guān)系 Fig. 1 Relationship between pitch angles and power

文本提出的風(fēng)電機(jī)組最優(yōu)槳距角辨識方法流程如圖2所示，主要分為4個部分：異常數(shù)據(jù)清洗、三維特性曲面擬合、最優(yōu)槳距角辨識和修正后發(fā)電量提升驗證。

圖2 風(fēng)電機(jī)組最優(yōu)槳距角辨識流程 Fig. 2 Identification process of optimal pitch angle of wind turbines

（1）異常數(shù)據(jù)清洗：為了消除異常數(shù)據(jù)點對整體分布規(guī)律和各變量間對應(yīng)關(guān)系的影響，需要對初始數(shù)據(jù)進(jìn)行預(yù)處理。首先剔除機(jī)組停機(jī)和故障時間段內(nèi)無意義的數(shù)據(jù)；然后基于改進(jìn)的DBSCAN聚類算法對“風(fēng)速-功率”“風(fēng)速-槳距角”等變量的離群數(shù)據(jù)及邊緣數(shù)據(jù)進(jìn)行識別并剔除；最后針對機(jī)組輸出功率的多影響因素耦合問題，對數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理，消除機(jī)型、量綱、空氣密度等因素的影響。

（2）三維特性曲面擬合：為了減小離散采樣數(shù)據(jù)誤差，首先提取預(yù)處理過后的機(jī)組“風(fēng)速-功率-槳距角”散點圖；然后，基于最小二乘法曲線曲面擬合方法，建立機(jī)組實際運行的“風(fēng)速-功率-槳距角”連續(xù)特性曲面。

（3）最優(yōu)槳距角辨識：為了辨識不同平均風(fēng)速下最優(yōu)槳距角，首先，通過LM（Levenberg- Marquardt）算法求解三維特性曲面變量函數(shù)關(guān)系；然后，以1 m/s為步長，通過辨識得到不同平均風(fēng)速下風(fēng)電機(jī)組最優(yōu)槳距角。

（4）最優(yōu)槳距角修正后的發(fā)電量提升驗證：為了驗證所提方法的有效性，首先在Bladed中對同一機(jī)組槳距角數(shù)值進(jìn)行修正并驗證；然后，將該機(jī)組模型提取到MATLAB中，在相同工況下再次進(jìn)行驗證。

2 基于運行數(shù)據(jù)的辨識算法

2.1 基于DBSCAN算法的離群數(shù)據(jù)剔除

考慮風(fēng)電機(jī)組主要特征量的相關(guān)性，選擇風(fēng)速、功率、槳距角3個變量進(jìn)行最優(yōu)槳距角辨識?！帮L(fēng)速-功率-槳距角”三維特性散點圖如圖3所示。

圖3 初始數(shù)據(jù)三維特性散點圖 Fig. 3 3-D characteristic scatter diagram of initial data

從圖3可以看出，運行數(shù)據(jù)中存在一定量的異常數(shù)據(jù)，其中包括離群數(shù)據(jù)及邊緣數(shù)據(jù)。異常數(shù)據(jù)產(chǎn)生的主要原因包括[7]：

（1）陣風(fēng)或者控制器延遲。此類異常數(shù)據(jù)量大，為非穩(wěn)態(tài)數(shù)據(jù)點，一般分布在正常數(shù)據(jù)簇周邊。

（2）風(fēng)電機(jī)組故障。此類異常數(shù)據(jù)較多，一般偏離正常數(shù)據(jù)且呈聚集分布。

（3）電磁干擾或者通信故障。此類異常數(shù)據(jù)量少，一般表現(xiàn)為離群數(shù)據(jù)。

在進(jìn)行機(jī)組特性曲面擬合之前必須清洗初始數(shù)據(jù)，否則異常數(shù)據(jù)的存在將會嚴(yán)重影響運行數(shù)據(jù)的整體分布，導(dǎo)致變量間函數(shù)關(guān)系誤差較大甚至錯誤。為了得到更準(zhǔn)確的數(shù)據(jù)，初始數(shù)據(jù)點中的邊界點和噪聲點也被看作是異常數(shù)據(jù)。

DBSCAN算法基本流程圖如圖4所示，算法中半徑Ep和鄰域密度閾值Mp是給定的。

圖4 基本DBSCAN算法流程圖 Fig. 4 Basic DBSCAN algorithm flow chart

本文對基本的聚類算法進(jìn)行改進(jìn)：結(jié)合初始數(shù)據(jù)特性，對上述2個重要參數(shù)進(jìn)行自動選擇，以確保聚類結(jié)果的正確性。

改進(jìn)的DBSCAN聚類算法主要的處理步驟為：

（1）根據(jù)初始數(shù)據(jù)量的大小，規(guī)定噪聲點的可接受占比σ。假定σ=2%；初始變化量k=2。

（2）令Mp=k。計算每個數(shù)據(jù)和其第k個最近對象間的距離d。初始數(shù)據(jù)中，所有的對象集合記為Dk。將Dk中累計概率位于6%～94%的數(shù)據(jù)記 為新的集合，以的數(shù)學(xué)期望作為Epk：

（3）統(tǒng)計初始數(shù)據(jù)里所有元素Epk鄰域內(nèi)點的個數(shù)，其數(shù)目集合記為Pk。同上，將累計概率 在6%～94%范圍內(nèi)的元素記為集合，將的數(shù) 學(xué)期望重新賦值給Mpk。

（4）以對初始數(shù)據(jù)進(jìn)行處理，并由聚類結(jié)果計算噪聲數(shù)據(jù)占比Rk。

（5）若k＞2，且|Rk-1-Rk|≤σ，則Mp=Mpk，Ep=Epk，流程結(jié)束；否則，令k=k+1，重復(fù)上述步驟。

選取“風(fēng)速-功率”“風(fēng)速-槳距角”等變量進(jìn)行聚類，結(jié)果如圖5所示。

圖5 基于改進(jìn)DBSCAN算法聚類結(jié)果 Fig. 5 Clustering results based on improved DBSCAN algorithm

除此之外，風(fēng)電機(jī)組性能在一定程度上會受到空氣密度、風(fēng)速等因素的干擾，因此需要考慮減小或消除此類因素的影響。

空氣密度和風(fēng)速標(biāo)準(zhǔn)化處理過程如下：

式中：ρ為空氣密度；B為氣壓；R為氣體常數(shù)，設(shè)置為287 J/(kg·K)；T為環(huán)境溫度；vn為標(biāo)準(zhǔn)空 氣密度下的風(fēng)速；v為實際風(fēng)速；0ρ為標(biāo)準(zhǔn)大氣 壓下的空氣密度，取1.225 kg/m3。

經(jīng)過數(shù)據(jù)預(yù)處理后的“風(fēng)速-功率-槳距角”三維特性散點圖如圖6所示。

圖6 預(yù)處理后三維特性散點圖 Fig. 6 Three dimensional characteristic scatter diagram after pretreatment

通過對比圖3、圖6可以看出，本文所采用改進(jìn)的DBSCAN算法可以有效剔除異常數(shù)據(jù)。清洗掉的異常值占45.3%。

2.2 基于最小二乘法機(jī)組特性曲面擬合

為了減小離散采樣數(shù)據(jù)誤差，需要對三維特性散點圖進(jìn)行擬合處理。

最小二乘法是一種廣泛被用來解決曲線曲面擬合問題的方法[8]。使用該方法可避免相鄰分段上的特性曲線不連續(xù)不平滑的問題，非常適合于數(shù)據(jù)量大且形狀復(fù)雜的離散數(shù)據(jù)擬合[9]。

最小二乘法擬合函數(shù)多項式可表示為：

式中：P為風(fēng)電機(jī)組輸出功率；v為風(fēng)速；β為槳距角；aj即為擬合曲面的一組系數(shù)。

為了得到不同平均風(fēng)速下風(fēng)電機(jī)組最優(yōu)槳距角，需要求解上述擬合函數(shù)。通常求解非線性方程的方法有梯度法和牛頓法。

梯度法的缺點是收斂速度慢，其迭代方程為：

式中：α為迭代步長；-gk表示負(fù)梯度方向。

牛頓法又稱海森矩陣法，其迭代方程為：

式中：Hk為海森矩陣。

與梯度法相比，牛頓法保留了泰勒級數(shù)的一階項且利用了二階項、考慮了梯度變化趨勢，其優(yōu)點是搜索方向更好、收斂速度快，缺點是因需要計算海森矩陣及其逆矩陣，所以其計算量大。

用雅可比矩陣近似代替海森矩陣，高斯-牛頓法很好地解決了這一問題。

考慮當(dāng)海森矩陣為不滿秩矩陣時則無法進(jìn)行迭代，故需引入一個單位矩陣：

將式（9）代入式（7），從而得到LM算法的迭代格式：

式中：μ表示阻尼因子，為正值。

LM算法解決了海森矩陣不滿秩無法迭代的問題，且具有梯度法和牛頓法的優(yōu)點。故本文采用LM算法求解擬合函數(shù)。

通過LM算法求解出在問題域上全部節(jié)點的系數(shù)，即可得到三維特性曲面擬合關(guān)系式：

“風(fēng)速-功率-槳距角”三維特性曲面擬合結(jié)果如圖7所示。

對于新式風(fēng)電機(jī)組，最優(yōu)槳距角的角度除0°外，還可以調(diào)節(jié)為正角度或負(fù)角度。本文以1 m/s為步長，在不同平均風(fēng)速下對基于最小二乘法擬合的曲面進(jìn)行切割，對切割形成的曲線求功率最大值對應(yīng)的槳距角。最終得到的不同平均風(fēng)速下風(fēng)電機(jī)組最優(yōu)槳距角如表1所示。

表1 不同平均風(fēng)速下風(fēng)電機(jī)組最優(yōu)槳距角 Tab. 1 Optimal pitch angle of wind turbine under different average wind speeds

3 算例驗證

本文以某5 MW風(fēng)電機(jī)組為對象在Bladed進(jìn)行建模。機(jī)組主要技術(shù)參數(shù)如表2所示。

表2 仿真機(jī)組主要技術(shù)參數(shù) Tab. 2 Main technical parameters of 5MW unit

仿真計算中，為模擬實際風(fēng)況，湍流風(fēng)平均風(fēng)速設(shè)置為4 m/s、6 m/s、8 m/s、10 m/s、12m/s。

為了模擬實際風(fēng)電機(jī)組最優(yōu)槳距角誤差情況，設(shè)置：最優(yōu)槳距角最大值為5°，最小值為-5°；每0.5°為1種工況條件，其中風(fēng)電機(jī)組理論槳距角為0.5°，故最優(yōu)槳距角工況共21種；每個工況時長600 s；每1 s輸出一條數(shù)據(jù)，共輸出數(shù)據(jù)61 950條數(shù)據(jù)。

輸出量包括風(fēng)速、功率、槳距角、發(fā)電機(jī)轉(zhuǎn)速以及轉(zhuǎn)矩。

仿真中出現(xiàn)的機(jī)組運行異?；蛱崆敖K止均不做處理，以此驗證本文所提方法的實用性。

首先，將槳距角數(shù)值設(shè)置為理論槳距角，仿真運行得到風(fēng)電機(jī)組輸出功率；然后計算得到機(jī)組發(fā)電量，視為“初始發(fā)電量”；同理，最優(yōu)槳距角仿真運行得到的發(fā)電量視為“修正后發(fā)電量”。

槳距角修正前后風(fēng)電機(jī)組1 h的發(fā)電量總體情況如圖8所示，具體數(shù)據(jù)如表3所示。

圖8 Bladed中槳距角修正前后發(fā)電量對比圖 Fig. 8 Comparison of power generation before and after pitch angle correction in bladed

表3 Bladed中槳距角修正前后發(fā)電量具體數(shù)據(jù) Tab. 3 Data of power generation before and after pitch angle correction in bladed

由表3及圖8可知，當(dāng)平均風(fēng)速在額定風(fēng)速以下時，通過本文提出的辨識方法求解得到的最優(yōu)槳距角能夠在一定程度上提升風(fēng)電機(jī)組發(fā)電量。平均風(fēng)速為7 m/s時，發(fā)電量提升最為明顯——相較于修改前提升了1.60%。當(dāng)平均風(fēng)速臨近額定風(fēng)速時，發(fā)電量的提升量較低。

將機(jī)組模型提取到MATLAB中，在MATLAB環(huán)境下進(jìn)一步驗證該最優(yōu)槳距角的有效性。針對風(fēng)電機(jī)組的強(qiáng)非線性，為獲得其線性化數(shù)學(xué)模型，考慮所需要的輸入、輸出以及模態(tài)。采用泰勒級數(shù)在各風(fēng)速點進(jìn)行線性化展開，可以獲得如式（12）所示風(fēng)電機(jī)組多輸入多輸出狀態(tài)空間方程。

式中：A、B、C、D表示三維系數(shù)矩陣；x表示狀態(tài)變量；u表示輸入變量；y表示輸出變量。

提取到模型后，在相同的工況下再次進(jìn)行發(fā)電量提升驗證。槳距角修正前后機(jī)組1 h的發(fā)電量情況如圖9所示，具體數(shù)據(jù)如表4所示。

圖9 MATLAB中槳距角修正前后發(fā)電量對比圖 Fig. 9 Comparison of power generation before and after pitch angle correction in MATLAB

表4 MATLAB中槳距角修正前后發(fā)電量具體數(shù)據(jù) Tab. 4 Data of power generation before and after pitch angle correction in MATLAB

由表4及圖9可知，在MATLAB仿真實驗中，辨識得到的最優(yōu)槳距角能夠提升風(fēng)電機(jī)組發(fā)電量。在平均風(fēng)速7 m/s的湍流風(fēng)工況下，發(fā)電量提升最為明顯——相較于修改前提升了1.59%。當(dāng)平均風(fēng)速臨近額定風(fēng)速時，發(fā)電量的提升量一般，平均風(fēng)速10 m/s的湍流風(fēng)工況下提升量為0.70%。

上述2種模型的計算結(jié)果均表明，在相同的工況下，最優(yōu)槳距角可提升風(fēng)電機(jī)組發(fā)電量，且提升量幾乎一致；這說明本文提出的最優(yōu)槳距角策略具有有效性。

4 結(jié)論

針對風(fēng)電機(jī)組最優(yōu)槳距角與理論槳距角不一致的問題，本文提出了一種基于運行數(shù)據(jù)驅(qū)動的最優(yōu)槳距角辨識方法。對同一參數(shù)的機(jī)組，將辨識得到的最優(yōu)槳距角分別通過Bladed和MATLAB進(jìn)行驗證。通過算例驗證分析可以獲得以下結(jié)論。

（1）針對算例驗證中5 MW風(fēng)電機(jī)組，通過本文方法可獲得不同平均風(fēng)速下最優(yōu)槳距角。當(dāng)平均風(fēng)速在額定風(fēng)速以下時，槳距角修正后能夠在一定程度上提升風(fēng)電機(jī)組發(fā)電量。在平均風(fēng)速7 m/s的湍流風(fēng)工況下，發(fā)電量提升最為明顯：相較于修改前提升1.59%～1.60%。

（2）在Bladed和MATLAB模型驗證中，槳距角修正后機(jī)組發(fā)電量的提升基本一致，故該辨識方法所獲得的最優(yōu)槳距角具有有效性。

本文所提出的風(fēng)電機(jī)組最優(yōu)槳距角辨識方法依賴于風(fēng)電機(jī)組運行數(shù)據(jù)。當(dāng)風(fēng)電機(jī)組模型參數(shù)發(fā)生改變時，本文結(jié)論中的數(shù)據(jù)不再適用。

基于此不足，進(jìn)一步研究方向為：考慮不同的風(fēng)電機(jī)組模型，驗證風(fēng)電機(jī)組最優(yōu)槳距角辨識方法的普遍性。