基于Makkonen結冰增長模型的風力機覆冰預報

李仲怡，葉庚姣，盧小鳳，李勇，唐百川

(1.廣西壯族自治區(qū)氣象災害防御技術中心，廣西 南寧 530022；2. 廣西壯族自治區(qū)氣象服務中心，廣西 南寧 530022；3.南寧海關緝私局，廣西 南寧 530022)

隨著風力發(fā)電技術日趨成熟，風能已成為具有開發(fā)潛力和競爭優(yōu)勢的新能源，近十年來得到了前所未有的關注和飛躍式發(fā)展。廣西風能資源較豐富，但是氣象災害頻發(fā)[1]。大部分風電場建設在高山丘陵地區(qū)，冬季氣溫較低、濕度大，極易發(fā)生覆冰災害且持續(xù)時間較長，風力機覆冰后會影響輸出功率，而且冰塊脫落可能造成安全事故[2-6]，嚴重的覆冰過程使風力機大面積停機，影響風電場經濟效益；因此，風電場對風力機覆冰起止時間和覆冰厚度的預報有迫切需求，以便于電場調度安排生產和提前采取防災措施。

由于發(fā)生載體不同，風力機覆冰與電線覆冰有本質區(qū)別。電線覆冰需超過承載負重才會影響輸電安全，而風力機一旦發(fā)生覆冰就會有影響，嚴重時甚至會持續(xù)一個月不能發(fā)電，因此需要區(qū)別對待，開展具有針對性的研究。國內外有大量關于電線覆冰的研究，其中Makkonen結冰增長模型是最典型且運用最為廣泛的模型之一[7]。Makkonen等于2000年提出電線覆冰模型，2001年將該模型應用于風力機覆冰模擬，研究了積冰模型對氣象變量的敏感性，并將模型預報結果與冰風洞覆冰試驗結果進行比較，結果應用于指導風力機葉片防冰加熱元件設計[8-11]。國內有大量運用Makkonen結冰增長模型開展電線覆冰的研究[12-14]，而關于風電場風力機覆冰的研究卻較少。由于缺乏有效的風電場風力機葉片覆冰厚度觀測數(shù)據，目前風力機覆冰的研究大多是基于結冰增長模型的仿真結果與冰風洞覆冰試驗結果來探討風力機結冰的影響因子[15-19]。受人工模型試驗條件限制，覆冰起止時間和厚度測量值存在誤差，較難應用于風電場風力機覆冰預報，因此開展實地風力機覆冰研究并建立本地化風力機覆冰預報模型對風電場有極大的應用價值。

結合實際風力機覆冰資料和歐洲中期天氣預報中心(European Centre for Medium-Range Weather Forecasts，ECMWF)模式預報資料，本文基于Makkonen結冰增長模型，首次以桂北風電場為研究對象對模型進行優(yōu)化，建立了廣西風力機覆冰預報模型，運用該模型對典型覆冰過程進行預報，并與實際風力機組覆冰起止時間比較，進行效果檢驗。

1 資料和方法

1.1 資料來源

數(shù)值模擬資料由ECMWF模式格點插值到風電場，ECMWF模式格點分辨率為12.5 km，時間分辨率為3 h，取前一日北京時間20時起報的未來72 h的預報數(shù)據驅動模型，預報要素包括降雨、氣溫、濕度、風速、氣壓等。氣象要素實況資料為風電場附近自動氣象觀測站數(shù)據，桂北2個風電場(風電場A和風電場B)提供風力機的覆冰起止時間、風力機參數(shù)等資料。

1.2 數(shù)據處理

2018年11月28日至2019年2月28日，桂北風電場A、B氣溫在0 ℃以下共有67個時次，其中ECMWF模式預報正確的只有47個時次。由于ECMWF模式對山區(qū)的氣溫預報誤差相對較大，通過與附近的自動站實況數(shù)據對比，采用線性回歸分析方法，得出實況氣溫和預報氣溫的關系，如式(1)；另外根據經驗，風電場A(B)氣溫預報值分別在1～3.2(3.8)℃、0～1.0 ℃時，實際氣溫已經分別在-1～0 ℃、-2～-1℃，因此對模式氣溫作式(2)的訂正。

{ti}℃=1.139+0.867 3{Ti}℃.

(1)

(2)

式(1)、(2)中：Ti為逐3 h氣溫預報值；ti為訂正以后的氣溫；風電場A的m取3.2 ℃，風電場B的m取3.8 ℃。

采用式(1)訂正ECMWF模式氣溫預報數(shù)據，氣溫平均絕對誤差降低了2.9%，0 ℃以下氣溫預報正確時次為30個，0 ℃以下氣溫平均絕對誤差反而增加了63.8%；采用式(2)訂正ECMWF模式氣溫預報數(shù)據，氣溫平均絕對誤差增加了9.6%，0 ℃以下氣溫預報正確時次為65個，0 ℃以下氣溫平均絕對誤差降低了23.9%。2019年1月風電場A的ECMWF模式氣溫訂正結果如圖1所示，采用經驗公式對0 ℃以下氣溫訂正效果較好的有221—224時次、282—283時次、312—317時次等。

圖1 2019年1月風電場A逐3 h氣溫曲線Fig.1 3 h temperature curves of wind farm A in January, 2019

總體來看，回歸分析方法對整體的氣溫訂正效果較好，但是對0 ℃以下氣溫訂正效果較差，而風力機覆冰預報模型對低溫預報準確率要求較高，由于經驗公式對0 ℃以下氣溫訂正效果較好，因此擇優(yōu)選擇式(2)對模式氣溫進行訂正處理。

1.3 結冰模型

由于本文側重研究積冰起止時間預報，不涉及冰形狀模擬，因此采用Makkonen結冰增長模型。覆冰過程主要由過冷卻水滴的碰撞、捕獲和凍結3個階段組成，積冰質量增長率表達式為

(3)

式中：M為積冰質量；α1為過冷卻水滴的碰撞系數(shù)；α2為捕獲系數(shù)；α3為凍結系數(shù)；w為粒子群含水量(質量分數(shù))；v為有效粒子速度(即粒子相對于風力機的速度)；A為水滴碰撞物體有效截面積；t為時間。

α1采用較常用的Finstad等[20]的計算方法：

其中，NFsys為電路的整體噪聲系數(shù)，NF1為第一級放大器的噪聲系數(shù)，G1為第一級放大器的增益，NF2為第二級放大器的噪聲系數(shù)。因為在實際的電路設計中，第一級放大器的設計對整體電路的噪聲系數(shù)影響最大，需要盡可能的降低第一級放大器的噪聲系數(shù)，同時提高第一級放大電路的增益，從而改善整體電路的噪聲特性。

(4)

式中：K=ρWd2v/9μD，φ=Re2/K；Re=ρadv/μ，K、φ量綱為一，Re為雷諾數(shù)；ρw為水密度，取1.0×103kg/m3；ρa為空氣密度，取1.293 kg/m3；v為相對風速；μ為空氣絕對黏度，取1.798 4×10-5kg/(ms)。D為電線直徑，Neil Davis等[21]以1 cm為槳葉厚度，模擬時假設直徑取0.144 cm圓柱體進行模擬；d為液滴直徑，F(xiàn)instad等[20]認為不用計算每個液滴直徑，只要知道中值體積直徑(median volume diameter，MVD)即可，周悅等[18]對其值dMV的估算方法為dMV=(3.672+μ1)/γ，其中普型參數(shù)μ1=min{15,1 000/Nc+2}，斜率﹝wWC為液態(tài)水含量(質量分數(shù))﹞，粒子總數(shù)濃度Nc取100 cm-3；Γ為伽馬函數(shù)。

α2采用Admirat等人提出的經驗公式：

(5)

對于凍結系數(shù)的計算，一般考慮干增長和濕增長。干增長時，所有被捕獲的液滴全部凍結在物體表面，α3=1；對于濕增長，被捕獲的液滴有部分沒有立即凍結，并且在積冰表面發(fā)生的熱傳遞過程中有部分積冰融化成為液態(tài)水，此時冰面熱平衡公式為

Qf+Qv=Qc+Qe+Ql+Qs.

(6)

式中：Qf為水滴凍結釋放的潛熱；Qv為氣流與冰面摩擦產生的熱；Qc為氣流帶走的熱感；Qe為冰面蒸發(fā)損失的熱；Ql為加熱過冷水滴到冰點損失的熱；Qs為短波輻射和長波輻射產生的熱。Makkonen等[8]通過解此方程，得到

(7)

式中：F=α1α2wv；E為表面未凍結部分(即液態(tài)水部分)占比，取0.3；h=Nu·λa/D，其中努塞爾數(shù)Nu采用謝真珍等[14]的計算方法，Nu=0.032Re0.85，空氣導熱系數(shù)λa=2.536 2×10-5J/(m·s·K)；σ為斯蒂芬波爾茲曼常數(shù)，取5.669 6×10-8W/(m2·K4)；a為輻射常數(shù)，取8.1×107K3；ε為水汽的摩爾分子比(即干空氣的分子質量與水蒸氣的分子質量的比值)，取0.62；p為氣壓，取實況值；ts、ta、td分別為冰面溫度、氣溫、液滴碰撞溫度，對于純水ts=0，對于云中積冰和凍雨結冰td=ta；水氣壓pa=p0kh/100，當T>0 ℃時{p0}Pa=611exp[17.27{T}℃/({T}℃+237.3)]，當T≤0 ℃時{p0}Pa=611exp[21.87{T}℃/({T}℃+265.5)]；k為電線表面局部恢復系數(shù)，取0.79；ps為飽和水氣壓，取6.17 Pa；水凍結潛熱Lf、蒸發(fā)潛熱Le分別取3.34×105J/(kg·K)、2.501×106J/(kg·K)；空氣比熱Cp、水的比熱Cw分別取1 004.07 J/(kg·K)、4 218 J/(kg·K)。

相對風速v采用Virk等[18]的計算方法，液態(tài)水含量wWC采用文獻[12]的計算方法，表達式為：

(8)

(9)

2 模型優(yōu)化

2.1 相對風速公式優(yōu)化

表1 GW82/1500風力機參數(shù)Tab.1 Parameters of GW82/1500 wind tubine

(10)

式中vf,i為逐3 h自然風速，i表示采樣點。當氣溫大于0 ℃時，不同自然風速下對應槳葉葉尖水汽粒子的相對風速如圖2所示。

圖2 氣溫大于0 ℃時相對風速與自然風速的關系曲線Fig.2 Relationship curve between relative wind speed and natural wind speed when the temperature is greater than 0 ℃

2.2 建立融冰影響因子

太陽輻射是影響風力機融冰的重要因素，但是ECMWF數(shù)值預報沒有這個要素，且風電場附近也缺乏太陽輻射的實況數(shù)據，而液態(tài)水含量能間接影響融冰，因此本文通過控制液態(tài)水含量來建立融冰影響因子。

當有降雨時，液態(tài)水含量wWC采用文獻[12]的計算方法，根據經驗當數(shù)值預報累計12 h無降水，并且預報氣溫大于0時，有利于融冰，液態(tài)水含量為0；當不能滿足較長時間無降水或氣溫低于0 ℃時，考慮空氣仍較濕潤，相對不利于融冰，液態(tài)水含量取0.1。表達式為：

(11)

式中TPi為降水預報值。

2.3 敏感性分析

影響風力機覆冰的外部因素是風力機的參數(shù)(見表1)，而內在影響因素有降水、相對濕度、氣溫、風速、氣壓、液態(tài)水含量等。為了驗證該模型在桂北風電機組的適用性，下面將針對氣溫、降水、風速進行敏感性分析。敏感性分析方案見表2，為方便分析把最大積冰厚度設置為100 mm。

表2 敏感性分析方案Tab.2 Sensitivity analysis schemes

氣溫、降水、風速影響覆冰厚度變化的曲線如圖3所示。圖3(a)所示為方案1，氣溫越低覆冰增長越快；隨著氣溫降低，相同的降溫幅度下，覆冰增長率的變化減小。圖3(b)所示為方案2，可見降水量越大，即液態(tài)水含量越高，覆冰速率越快。圖3(c)所示為方案3，隨著風速的增大，覆冰增長反而越慢，因為相對風速的增加會增大碰撞系數(shù)，但是翼型表面的熱力學平衡過程有2個對立影響：一是風速增大物體表面的對流散熱，利于覆冰量增加；二是風速增大，物體表面摩擦生熱以及水滴動能增大，利于覆冰量減少[16]。因此當風速在適宜范圍內時，覆冰增長速度才較大。

圖3 不同條件下覆冰厚度變化曲線Fig.3 Variation curves of icing thickness of wind turbine under different conditions

3 2019年冬季覆冰預報檢驗

將ECWMF氣象預報模式的降水、相對濕度、氣溫、風速、氣壓等資料代入建立的風力機覆冰預報模型，對桂北風電場A、B2019年11月1日至2020年3月31日的覆冰過程進行逐日預報，其中風電場A、B分別有66臺、340臺風力機，預報檢驗結果見表3、表4。覆冰停機時間指的是由覆冰引起的第1臺風力機停機的時間；融冰開機時間指風力機全部融冰并正常運行的時間；覆冰過程報出率是指風電場風力機覆冰起止時間不漏報的概率。

由表3可見，風電場A共8次過程，提前24 h、48 h、72 h覆冰過程報出率分別為100%、100%、75%。提前24 h、48 h、72 h覆冰開始時間預報的平均絕對誤差分別是3.8 h、4.1 h、3.7 h，預報誤差小于5 h的過程占比分別為75%、75%、83%，預報誤差小于2 h的過程占比分別為38%、25%、17%。風電場A融冰開機時間受金中旺線輸電線路融冰影響，第5次過程的覆冰結束預報時間誤差較大，因此剔除異常偏大過程，提前24 h、48 h、72 h預報的覆冰結束時間平均絕對誤差分別是7.8 h、7.8 h、13.2 h，預報誤差小于5 h的過程占比分別為29%、43%、14%，預報誤差小于2 h的過程占比分別為29%、43%、14%。

表3 風電場A風力機覆冰預報檢驗結果Tab.3 Test results of wind turbine icing prediction in wind farm A

由表4可見，風電場B共9次過程，提前24 h、48 h、72 h覆冰過程報出率分別為100%、100%、67%。提前24 h、48 h、72 h預報的覆冰開始時間平均絕對誤差分別是4.8 h、4.7 h、5.6 h，預報誤差小于5 h的過程占比分別為67%、78%、67%，預報誤差小于2 h的過程占比分別為33%、44%、50%。風電場B提前24 h、48 h、72 h預報的覆冰結束時間平均絕對誤差分別是6.1 h、11.6 h、11.7 h，預報誤差小于5 h的過程占比分別為67%、44%、43%，預報誤差小于2 h的過程占比分別為44%、11%、29%。

表4 風電場B風力機覆冰預報檢驗結果Tab.4 Test results of wind turbine icing prediction in wind farm B

總體來看，風電場預報覆冰開始時間早于實際停機時間，融冰開機時間反之，因為風力機可能達到一定覆冰厚度后才會出現(xiàn)輸出功率異常或直接停機的情況，據統(tǒng)計覆冰停機時覆冰厚度為0.5～1 mm。預報的覆冰厚度越大，覆冰時間越長，風電機組停機臺數(shù)越多；因此，可以根據模型預報的覆冰厚度和時長來預判覆冰過程強度，預估停機數(shù)量。另外，風力機覆冰厚度較小時，可能恢復正常運轉，據統(tǒng)計融冰開機時覆冰厚度為0.5～2 mm。如表3中序號1過程，風力機停機時間僅4.2 h，風力機僅停機2臺，這種覆冰時間短、過程強度較小的風力機覆冰過程誤差相對較大，原因可能是ECWMF模式對較弱降溫過程的氣溫預報誤差較大。

從以上預報結果可以看出提前48～72 h的風力機覆冰過程預報結果與實際覆冰情況相近，預報風力機覆冰起止時間大致相符，能為風電場合理安排生產和采取防災措施提供有效的參考。

4 結論

本文基于Makkonen結冰增長模型，通過優(yōu)化模型和氣象因子敏感性分析，建立了風力機覆冰預報模型，采用優(yōu)化的風力機覆冰預報模型對桂北風電場A、B 2019年11月1日至2020年3月31日期間的覆冰過程進行了預報和檢驗，結論如下：

a)通過訂正氣象數(shù)值預報的氣溫值、優(yōu)化相對風速計算方法和建立融冰影響因子，使風力機覆冰預報模型更符合風力機運行和覆冰特點，以提高適用性和預報準確率。

b)敏感性分析結果表明，在適宜的風速下，溫度越低，降水越大，覆冰增長越快，符合覆冰物理變化規(guī)律。

c)該模型對風電機組覆冰時間預報具備預報時效長、預報準確率高的特點。48 h時效內覆冰過程報出率為100%，覆冰起止時間預報誤差分別為3～5 h、6～12 h；72 h時效內風力機覆冰預報效果較穩(wěn)定，65%以上的覆冰開始時間預報誤差可控制在5 h內。風電場可以根據預報的覆冰起止時間、覆冰時長和覆冰厚度等預判覆冰過程強度，提前48～72 h安排生產和提前啟動防災預案。

風力機結冰對風電場經濟和安全造成嚴重影響，本文基于Makkonen結冰增長模型，建立了風力機覆冰預報模型，研究結果可對整個風電場或者分區(qū)結冰提供預報參考。由于缺乏風電場風力機覆冰厚度觀測數(shù)據，目前風力機覆冰模型的厚度預報無法檢驗。本文模型預報的過程最大覆冰厚度大小與覆冰造成的停機數(shù)量有較好的對應關系，在今后的研究中，可以根據覆冰過程停機情況和覆冰厚度建立關系，表征覆冰過程強度，或者在風電場搭建風力機葉片覆冰觀測設備。驅動模式ECWMF模式分辨率為12.5 km，針對單臺風力機的預報，還需要引進更高分辨率的預報模式，以提高預報精細化程度。