基于威布爾分布的風功率密度計算方法比較

李化

（潤陽能源技術有限公司,天津 300300）

0 引言

風資源評估是風電開發的重要環節[1]，風功率密度又是風能資源評估的重要參數之一[2]，是決定風電場是否具有開發潛能的重要評判標準。在年平均風速和空氣密度相同的情況下，不同的風頻分布將計算得出不同的風功率密度，因而準確分析風資源的風頻分布具有重要意義。

目前，對風頻分布的研究多采用威布爾分布，它能比較好地擬合風資源的實際情況。劉鵬等[3]采用了矩估計法模擬風頻威布爾分布，此方法簡單方便，缺點是不能完全利用風數據樣本信息，模擬精度不高。Justus 和Lysen 提出的經驗計算方法[4]，被2002版《全國風能評價技術規定》采納，因簡單、實用、耗時短，在風電場建設初期做初步評估時被廣泛采用。Akdag 和Dinler 提出了用能量因子法來模擬威布爾分布參數[5]，此方法目前被廣泛應用于風資源流體力學評估模型軟件中，如Meteodyn Universe。Stevens MJM 將最大似然法運用于風資源風頻分布的模擬計算，并以此估算風電場發電量[6]。最小二乘法由法國勒讓德創立，是處理實測值的純粹代數方法，此方法被廣泛應用于各行業，這是一種尋找最佳解，揭示最接近真實情形的狀態[7]。

估算風功率密度的方法較多，準確地計算風功率密度有賴于風頻威布爾分布擬合的準確性[8]。但是，目前還缺少基于威布爾分布擬合準確性計算風功率密度方面的研究。因此，文章對目前國際常用的5 種風頻威布爾擬合方法進行對比研究，采用決定系數對威布爾擬合結果準確度進行評價，并在此基礎上利用威布爾分布計算風功率密度，并與實測風數據計算所得風功率密度相比較絕對誤差和相對誤差，從而確定哪些方法更可靠，更準確。文章選取我國不同經緯度和地形地貌的4 座測風塔實測數據進行威布爾分布擬合和決定系數分析，并計算相應的風功率密度，用案例進行對比驗證得出較優方法，可為工程師在實際的風資源評估工作中提供參考和指導。

1 威布爾分布和風功率密度

威布爾分布是由瑞典科學家威布爾從材料強度的統計理論中推導出來的，是日常生活中一種常用的分布，在風電領域的風頻分布中也被廣泛采納應用。

威布爾的形狀參數 k值反映了風速分布的范圍寬窄，一般情況下，k值越小，風速分布越寬。當0

式中：

F(v)——威布爾風頻累計概率密度（%）；

v ——風速（m/s）；

c ——威布爾尺度參數（m/s）；

k ——為威布爾形狀參數。

對威布爾累計概率密度函數（1）求導數，得到威布爾風頻分布概率密度函數，如式（2）所示：

式中：

f(v)——威布爾風頻分布概率密度（%）。

風功率密度是指空氣流垂直通過單位截面積的風能，是表征1 個地方風能資源多少的1 個指標，在風資源評估中一般指年均風功率密度，某地的風功率密度越大，說明該地區風能資源越好[9]。風功率密度計算方法有兩種，一種是基于時間序列實測風數據[10]，風功率密度與風速的3 次方成正比，其計算公式如式（3）所示：

式中：

ρ——空氣密度（kg/m3）；

vi——風速時間序列值（m/s)；

n——時間序列。

另一種是選擇合適風速分布的數學模型來模擬風速變化，常見的有雙參數威布爾分布、瑞利分布、伽馬分布、對數正態分布，大多學者認為雙參數威布爾分布可以較好的模擬風速變化[11-13]?；谕紶柗植嫉娘L功率密度計算公式如式（4）所示：

2 威布爾分布參數擬合方法及精度驗證

2.1 EMJ Justus 經驗方法

該經驗方法由Justus 提出[5]，我國頒布的《全國風能評價技術規定》中也采用此種方法，此方法計算簡單。其 k 值和 c值計算公式如式（5）和（6）[14]所示：

式中：

σ ——風速標準差（m/s）；Γ（1+1/k）——伽馬函數；

2.2 EML Lysen 經驗方法

此方法由Lysen 提出[15]，在Justus 提出的方法上有所改進，其 k 值計算由式（5）給出，c值計算公式如式（7）所示：

2.3 EPF 能量因子方法

基于能量守恒規律，使其威布爾擬合計算得到的風功率密度與實測數據計算所得的風功率密度差值盡可能小的情況下求 k 和 c值。此方法首先求出能量因子 Epf，其計算公式如式（8）所示，再計算 k 值[16]，k值計算公式如式（9）所示，c值計算公式同式（7）：

式中：

Epf——風功率密度能量因子（W/m2）。

2.4 MLE 最大似然法

1821 年由德國數學家高斯提出最大似然法，此方法首次用于計算發電量可以追溯到Stevens MJM[6]，用這種模型來估算風頻分布參數 k 和 c值，其中蘊含的規則是：在求得的 k 和 c這個特定的威布爾分布函數中，能得到現在的實測風數據序列樣本概率最大。

假設某一風速時間序列樣本為：

威布爾概率密度似然函數 f 在某一風速 vi下的公式如式（10）所示：

則風速時間序列樣本 V 產生的聯合密度函數為：

為了計算方便，等式兩邊取對數：

上述方程分別對k，c求導，并令其分別為0，見式（13）和式（15），即可求得威布爾分布的參數 k 和c[17]，具體見式（14）和式（16）：

2.5 LLSA 最小二乘法

風速時間序列某一風速的累計概率密度分布函數公式如式（17）所示：

為計算方便，等式兩邊取對數如式（18）所示：

累計概率密度函數經過上述的變化推導就變成常見的線性 y=mx+b形式。最小二乘法的理論基石是實測值與模擬求得的結果值之間的差的平方總和最小。在用于威布爾擬合中，其公式如式（19）所示：

式中：

SSerr——為實測值與模擬結果值差的平方總和。

式（19）對 k 和 c 分別求導，并分別令其為0,得出k值和 c值[18]，如式（20）和式（21）所示：

上述最大似然法和最小二乘法所涉及的公式屬于理論推導，運算較復雜，在涉及數據繁多的測風數據處理時，常常借助Excel 中的規劃求解器Slover 來解決，Slover 規劃求解器是一種多次運算組合模擬參數使得目標函數取得最大值或者最小值的解決辦法[19],此方法方便快捷，批量處理測風數據能提高計算的效率和準確性。具體運用于此文章時將多次組合模擬 k 和 c值，使得聯合密度函數取得最大值，最小二乘法函數即式（19）中的 SSrr取得最小值，此時相應的 k 和 c即為最優解。

2.6 威布爾模擬精度驗證

在實際運用中，針對風場實測數據的風頻分布情況一般采用上述5 種方法進行擬合，文章采用決定系數R2來衡量其擬合結果的精確度[20]，其計算公式如式（22）所示：

式中：

fw(vi) ——威布爾分布參數計算風速為 vi時的風頻（%）；

fm(vi) ——實測風速序列為 vi時的風頻（%）；——實測風速序列風頻的平均值（%）；

0.9

0.8

0.7

R2≤0.7，擬合精度差。

3 案例分析

文章選取我國不同經緯度，不同地形地貌及風資源條件各不相同的4 座測風塔，1 個完整年統計值為10 min 的風速樣本52 560 個數據進行分析[21]。測風塔1 位于內蒙古茂明的平坦地形，植被稀少地貌非常簡單；測風塔2 位于云南昭通東部海拔為2 100 m 左右的多山地帶，植被較多，地形復雜；測風塔3位于新疆小草湖的半沙漠干旱地帶；測風塔4 位于廣東低緯度甲東地區的沿海灘涂附近。其測風塔所在地理區域位置如圖1 所示。

圖1 4 座測風塔地理位置圖Fig.1 The geographic position of 4 masts

測風塔的經緯度，所在地空氣密度，年平均風速和實測風功率密度等基本參數如表1 所示。

將測風數據以0.5 m/s 作為一個風速區間范圍進行頻率統計，采用上文介紹的5 種威布爾擬合方法分別對收集到的4 座測風塔數據進行模擬計算得k 和 c值，用威布爾分布函數可以求出各風速段相對應的頻率，在此基礎上利用式（4）計算威布爾分布風功率密度和式（22）計算決定系數，各測風塔結果分別如表2～表5 所示。

表2 測風塔1 威布爾分布和模擬風功率密度Tab.2 Simulation results of Weibull distribution and power density of mast 1

表3 測風塔2 威布爾分布和模擬風功率密度Tab.3 Simulation results of Weibull distribution and power density of mast 2

表4 測風塔3 威布爾分布和模擬風功率密度Tab.4 Simulation results of Weibull distribution and power density of mast 3

表5 測風塔4 威布爾分布和模擬風功率密度Tab.5 Simulation results of Weibull distribution and power density of mast 4

從表2～表5 模擬結果可知，4 座測風塔用經驗法EMJ 和EML 計算的參數值K 分別相等，C 值均相差0.01，非常接近。

測風塔1 和測風塔4 用5 種威布爾方法模擬風頻結果均好，決定系數均高于0.90，能量因子法EPF和最大似然法MLE 的決定系數相對更高，均高于0.96。相應威布爾計算所得的風功率密度與實測數據計算所得風功率密度絕對誤差和相對誤差也小，其中測風塔1 使用能量因子法EPF 和最大似然法MLE 與實測風數據計算所得風功率467.50 W/m2絕對誤差分別為9.70 W/m2和8.00 W/m2，相對誤差值分別為2.1%和1.7%。測風塔4 用5 種威布爾方法計算所得風功率密度和實測風功率密度絕對誤差分別為46.12 W/m2，44.44 W/m2，5.89 W/m2，3.66 W/m2和66.88 W/m2。

測風塔2 和測風塔3 威布爾模擬的決定系數橫向比較也有相似的趨勢，能量因子法EPF 和最大似然法MLE 的決定系數高于其他3 種方法，其計算所得的風功率密度與實測數據所得風功率密度更加接近，絕對值誤差和相對誤差均小于其他3 種方法。盡管測風塔3 模擬決定系數整體偏低，其分布范圍為0.69～0.85 之間，但能量因子法EPF 和最大似然法MLE 結果擬合精度較好分別為0.83 和0.85，其相對較高的趨勢仍未變化。

測風塔3 用能量因子法EPF 和最大似然法MLE 計算所得的風功率密度與實測數據計算所得風功率密度絕對誤差分別為28.81 W/m2和16.53 W/m2，相對誤差分別為4.7% 和2.7%；經驗法EMJ和EML 絕對誤差較大，分別為125.94 W/m2和127.93 W/m2，最小二乘法LLSA 為68.46 W/m2，相對誤差分別為20.6%，17.3%和11.2%。

測風塔3 用經驗法EMJ 和EML 及最小二乘法LLSA 威布爾模擬結果一般，決定系數分別為0.69，0.69 和0.75，在此基礎上利用風頻分布計算所得的風功率密度誤差更大，這也更好地驗證了風頻分布擬合準確與否對風功率密度計算準確度影響至關重要，風速與風功率密度是3 次方關系。

4 座測風塔實測風頻分布圖和威布爾分布模擬分布情況如圖2～圖5 所示。

圖2 測風塔1 實測風速風頻和威布爾擬合分布圖Fig.2 Wind speed frequency of the measurement data and Weibull distribution simulation of mast 1

圖3 測風塔2 實測風速風頻和威布爾擬合分布圖Fig.3 Wind speed frequency of the measurement data and Weibull distribution simulation of mast 2

圖4 測風塔3 實測風速風頻和威布爾擬合分布圖Fig.4 Wind speed frequency of the measurement data and Weibull distribution simulation of mast 3

圖5 測風塔4 實測風速風頻和威布爾擬合分布圖Fig.5 Wind speed frequency of the measurement data and Weibull distribution simulation of mast 4

4 結論

文章針對風頻分布擬合方法進行比較研究和風功率密度計算結果驗證，引用決定系數對風頻分布擬合度進行優劣評價，確定威布爾擬合常用的5 種方法準確性高低。為了使結論具有普適性和代表性，文章采用我國不同經緯度，地形地貌各異和風資源稟賦程度各不相同的4 座實測測風塔1 年統計數據進行驗證。

最終得出結論：威布爾擬合常用的5 種計算方法中，能量因子法EPF 和最大似然法MLE 能更好地擬合實測風數據的風頻，采用威布爾計算所得的風功率密度和實測風數據計算的風功率密度更為接近，誤差較小，決定系數值較大，其風頻的擬合度優于經驗法和最小二乘法；經驗法和最小二乘法擬合度相對差一些,相應計算的風功率密度誤差也大一些。盡管某些測風塔用5 種方法擬合所得的決定系數相對都較小，擬合度一般，但橫向比較能量因子法EPF 和最大似然法MLE 的擬合度仍較經驗法和最小二乘法有優勢，相應計算所得的風功率密度誤差仍小于其他3 種方法。

這一結果與目前我國較流行的流體力學風資源模擬軟件Meteodyn Universe 的運算規則相吻合[22]，這款軟件在做各高層繪圖計算時用的是能量因子法EPF，在各風機機位點參數計算 k 值，c值時運用的是最大似然法MLE。

因此，在實際項目風資源評估中，建議采用能量因子法EPF 和最大似然法MLE 來模擬風數據的風頻分布，計算得出更加接近實際的風功率密度，可為更好地了解掌握風場的風資源提供可靠依據。