劉福國(guó) 國(guó)欽光，殷炳毅，王守恩

(1.國(guó)網(wǎng)山東省電力公司電力科學(xué)研究院，山東 濟(jì)南 250003；2.山東電力研究院，山東 濟(jì)南 250003；3.華電濰坊發(fā)電有限公司，山東 濰坊 261000)

1 引 言

檢測(cè)發(fā)電機(jī)組效率時(shí)，鍋爐燃煤發(fā)熱量還無(wú)法連續(xù)測(cè)量，只能通過(guò)采集入爐煤樣品，在實(shí)驗(yàn)室化驗(yàn)分析得到[1]。用少量樣品煤的化驗(yàn)數(shù)據(jù)表示檢測(cè)周期內(nèi)入爐煤發(fā)熱量，其不確定性來(lái)自樣品采集和實(shí)驗(yàn)室分析等2個(gè)環(huán)節(jié)[2]，采樣環(huán)節(jié)引起的不確定度占95%以上[3]。文獻(xiàn)[4]研究了鍋爐采樣樣品煤發(fā)熱量的概率分布特性，對(duì)單個(gè)樣品發(fā)熱量的不確定度進(jìn)行了評(píng)定。在實(shí)際生產(chǎn)中，機(jī)組效率常用供電煤耗表示，電廠日常供電煤耗檢測(cè)周期有時(shí)長(zhǎng)達(dá)1個(gè)月、1個(gè)季度甚至1年。檢測(cè)周期內(nèi)的燃煤發(fā)熱量一般采用多個(gè)樣品煤發(fā)熱量的平均值[5]。因此，多次采樣樣品發(fā)熱量均值的不確定度評(píng)定對(duì)機(jī)組效率檢測(cè)具有重要意義。

μx=μ

(1)

(2)

當(dāng)總體N為正態(tài)分布或樣本容量n>30時(shí)，樣本均值的抽樣分布可作為正態(tài)分布，正態(tài)分布的參數(shù)μx、σx根據(jù)式(1)、式(2)計(jì)算，故可對(duì)n個(gè)樣品發(fā)熱量均值的不確定度進(jìn)行評(píng)定。

但當(dāng)總體N為非正態(tài)分布且樣本容量n<30時(shí)，樣本均值的抽樣分布為非正態(tài)分布，此時(shí)，如何對(duì)n個(gè)樣品發(fā)熱量均值的不確定度進(jìn)行評(píng)定，是本文要解決的問(wèn)題。

2 虛擬采樣及發(fā)熱量均值不確定度

非正態(tài)總體N通常來(lái)自實(shí)際采樣實(shí)驗(yàn)，為已知數(shù)據(jù)。要對(duì)n個(gè)樣品發(fā)熱量均值的不確定度進(jìn)行評(píng)定，首先獲取樣本均值的抽樣樣本，傳統(tǒng)方法[9]是從總體的N個(gè)發(fā)熱量中隨機(jī)抽取一個(gè)容量為n的樣本，計(jì)算該樣本均值，采用同樣方式進(jìn)行多次抽樣，得到樣本均值的抽樣樣本，利用該樣本對(duì)均值的不確定度進(jìn)行評(píng)定。從總體的N個(gè)發(fā)熱量中抽取容量為n的樣本，重復(fù)抽樣總計(jì)可得到Nn個(gè)樣本，因此樣本均值抽樣樣本的最大容量為Nn。

獲得樣本均值的抽樣樣本的另一種方法是基于蒙特卡洛法[10，11]的隨機(jī)采樣實(shí)驗(yàn)，與傳統(tǒng)抽樣方法相比，隨機(jī)采樣試驗(yàn)獲得的虛擬樣本容量是任意的，采用隨機(jī)采樣試驗(yàn)對(duì)n個(gè)樣品發(fā)熱量均值的不確定度進(jìn)行評(píng)定的步驟為：

(1) 計(jì)算總體N的發(fā)熱量概率密度。

(2) 將總體N發(fā)熱量概率密度表示成多項(xiàng)式函數(shù)。

(3) 建立虛擬采樣系統(tǒng)。按照總體N的概率密度產(chǎn)生隨機(jī)數(shù)，作為虛擬采樣樣品的發(fā)熱量。這種利用給定概率密度函數(shù)進(jìn)行隨機(jī)樣抽的系統(tǒng)，看作是樣品發(fā)熱量的虛擬采樣系統(tǒng)。

(4) 利用虛擬采樣系統(tǒng)進(jìn)行n次采樣，得到n個(gè)虛擬樣品的發(fā)熱量數(shù)據(jù)，計(jì)算發(fā)熱量的平均值，得到第1個(gè)發(fā)熱量均值；采用同樣的方法得到m個(gè)發(fā)熱量均值，作為發(fā)熱量均值的抽樣樣本。

(5) 計(jì)算發(fā)熱量均值樣本的平均值、標(biāo)準(zhǔn)差及概率密度分布，在給定置信概率下求取包含系數(shù)，得到n個(gè)樣品發(fā)熱量均值擴(kuò)展不確定度。

以某電廠鍋爐為例，對(duì)入爐煤發(fā)熱量均值不確定度進(jìn)行評(píng)定。

2.1 樣本N的獲得

表1為某電廠2017年5～7月入爐煤實(shí)際采樣試驗(yàn)樣本N，該樣本包含276個(gè)樣品煤的發(fā)熱量[12]。這些發(fā)熱量數(shù)據(jù)既包含實(shí)驗(yàn)室化驗(yàn)分析引起的不確定性，又包含采制樣本環(huán)節(jié)引起的不確定性。研究表明，實(shí)驗(yàn)室分析環(huán)節(jié)引起的不確定度約為0.05 MJ/kg[13～15]，而采制樣本環(huán)節(jié)引起的不確定度為0.5～2.5 MJ/kg[2，4，12]。因此，忽略實(shí)驗(yàn)室分析環(huán)節(jié)的不確定度，可認(rèn)為表1中發(fā)熱量數(shù)據(jù)的不確定性是由采制樣本引起的。利用表1的數(shù)據(jù)對(duì)n次采樣發(fā)熱量均值的不確定度進(jìn)行評(píng)定。

表1給出的樣本平均值μN(yùn)=19.686 MJ/kg，標(biāo)準(zhǔn)差σN=0.933 MJ/kg。

表1 總體樣本N的發(fā)熱量Tab.1 Calorific value of sample N MJ/kg

2.2 樣本N的概率密度函數(shù)

總樣本N中，發(fā)熱量最小值為16.894 MJ/kg，最大值為21.217 MJ/kg。將此區(qū)間等分成若干個(gè)子區(qū)間，統(tǒng)計(jì)每個(gè)子區(qū)間內(nèi)的樣品個(gè)數(shù)，計(jì)算落入?yún)^(qū)間的概率，概率密度是概率的變化率。利用Matlab統(tǒng)計(jì)工具箱中的hist函數(shù)可以統(tǒng)計(jì)落入每個(gè)子區(qū)間內(nèi)的樣品數(shù)，利用ksdensity函數(shù)可以對(duì)樣本N的概率密度進(jìn)行核心平滑估計(jì)[16]。由圖1給出的通過(guò)區(qū)間統(tǒng)計(jì)得到的概率密度的直方圖以及概率密度平滑估計(jì)可以看出，樣本發(fā)熱量不服從正態(tài)分布。概率密度平滑估計(jì)結(jié)果可用8階多項(xiàng)式函數(shù)f表示為：

圖1 總體樣本N的概率密度Fig.1 Probability density of population sample N

(3)

對(duì)于總體樣本N，式(3)所示的多項(xiàng)式函數(shù)由分段函數(shù)f1和f2組成，在不同區(qū)間上，函數(shù)f1和f2的系數(shù)見(jiàn)表2。

表2 多項(xiàng)式函數(shù)f的系數(shù)Tab.2 Coefficient of polynomial function f

2.3 虛擬采樣

根據(jù)概率密度函數(shù)式(3)，采用挑選抽樣方法產(chǎn)生發(fā)熱量樣本，這種方法由馮·諾依曼在1951年提出[17]。如圖2所示，當(dāng)已知概率密度函數(shù)f(Q)時(shí)，利用挑選抽樣法產(chǎn)生發(fā)熱量隨機(jī)數(shù)的主要步驟為：

(1) 在發(fā)熱量取值區(qū)間[a，b]內(nèi)，產(chǎn)生均勻隨機(jī)數(shù)Q。

Q=(b-a)r1+a

(4)

式中：r1為區(qū)間[0，1]均勻隨機(jī)數(shù)，r1∈[0，1]。

(2) 在[0，fmax]區(qū)間內(nèi)，產(chǎn)生均勻隨機(jī)數(shù)f。

f=r2fmax

(5)

式中：r2為區(qū)間[0，1]均勻隨機(jī)數(shù)，r2∈[0，1]；fmax為函數(shù)f(Q)在區(qū)間[a，b]上的最大值。

(3) 當(dāng)f≤f(Q)時(shí)，接受Q為所選的隨機(jī)數(shù)，否則，返回第(1)步，重新抽取1對(duì)(Q，f)。

根據(jù)圖2，挑選樣抽法的幾何解釋為：隨機(jī)選取位于矩形abcd內(nèi)的點(diǎn)(Q，f)，選擇位于曲線f(Q)下面的點(diǎn)，它們對(duì)應(yīng)的發(fā)熱量Q服從概率密度為f(Q)的分布。

圖2 根據(jù)給定概率密度函數(shù)f(Q)產(chǎn)生隨機(jī)數(shù)QFig.2 Generating random numbers Q based on a given probability density function f(Q)

2.4 均值樣本的獲得及評(píng)定

要對(duì)n次采樣發(fā)熱量均值的不確定度進(jìn)行評(píng)定，首先利用第2.3節(jié)的虛擬采樣技術(shù)獲得均值樣本。按式(3)給定的概率密度函數(shù)，產(chǎn)生n個(gè)發(fā)熱量隨機(jī)數(shù)，計(jì)算它們的平均值，作為發(fā)熱量均值的第1個(gè)采樣數(shù)據(jù)。用相同方式得到m個(gè)發(fā)熱量均值的采樣數(shù)據(jù)，作為發(fā)熱量均值樣本。

取n=4和n=50，采用虛擬采樣方法分別得到發(fā)熱量均值樣本，樣本容量m取為m=10 000。對(duì)樣本進(jìn)行統(tǒng)計(jì)計(jì)算，n=4時(shí)，樣本平均值μ4=19.689 MJ/kg，標(biāo)準(zhǔn)差σ4=0.480 MJ/kg；n=50時(shí)，樣本平均值μ50=19.689 MJ/kg，標(biāo)準(zhǔn)差σ50=0.135 MJ/kg。利用ksdensity函數(shù)對(duì)2個(gè)樣本的概率密度進(jìn)行核心平滑估計(jì)，n=4時(shí)，結(jié)果見(jiàn)圖3，其概率密度可擬合成式(3)所示的多項(xiàng)式函數(shù)f(Q)，限于篇幅，這里不再給出多項(xiàng)式系數(shù)；n=50時(shí)，區(qū)間統(tǒng)計(jì)得到的概率密度直方圖以及函數(shù)ksdensity的平滑估計(jì)結(jié)果見(jiàn)圖4，發(fā)熱量均值樣本非常接近正態(tài)分布，圖4中還給出了正態(tài)分布N(19.686，0.1352)的概率密度，可以看出，它和平滑估計(jì)結(jié)果吻合良好。

圖3 n=4時(shí)發(fā)熱量平均值的概率密度分布Fig.3 Probabilistic density distribution of the mean calorific value at n=4

圖4 n=50時(shí)發(fā)熱量平均值的概率密度分布與正態(tài)分布的對(duì)比Fig.4 Comparison of probability density of average calorific value with normal distribution at n=50

如圖3，對(duì)n=4時(shí)發(fā)熱量均值的不確定度進(jìn)行評(píng)定，就是求給定的置信概率ε下的包含因子k，根據(jù)置信概率的定義[18]得到關(guān)于包含因子k的方程如式(6)。

求解式(6)，可得到包含因子k。求解時(shí)取置信概率ε=0.95，式(6)中μ4=19.689，σ4=0.480。采用迭代法求解式(6)，得到k=1.934。因此，n=4時(shí)發(fā)熱量均值的擴(kuò)展不確定度為kσ4=0.928 MJ/kg，此時(shí)發(fā)熱量均值結(jié)果可表示為19.689±0.928 MJ/kg。

(6)

根據(jù)虛擬采樣方法，同樣可對(duì)n=50時(shí)發(fā)熱量均值的不確定度進(jìn)行評(píng)定。利用Matlab的jbtest函數(shù)對(duì)n=50時(shí)發(fā)熱量均值樣本進(jìn)行正態(tài)性檢驗(yàn)[19]。結(jié)果表明，均值樣本服從正態(tài)分布，置信概率ε=0.95對(duì)應(yīng)的包含因子k=1.96[20,21]，發(fā)熱量均值的擴(kuò)展不確定度為kσ50=0.265 MJ/kg。因此，n=50時(shí)發(fā)熱量均值的采樣化驗(yàn)結(jié)果可表示為19.689±0.265 MJ/kg。

3 結(jié) 論

(1) 在發(fā)電機(jī)組效率檢測(cè)時(shí)，鍋爐燃煤多次采樣樣品發(fā)熱量均值的不確定度評(píng)定對(duì)于機(jī)組效率檢測(cè)具有重要意義。

(2) 已知燃煤發(fā)熱量樣本，要對(duì)n個(gè)樣品發(fā)熱量均值的不確定度進(jìn)行評(píng)定，當(dāng)發(fā)熱量樣本正態(tài)分布或n>30時(shí)，則n個(gè)樣品發(fā)熱量均值服從正態(tài)分布，該正態(tài)分布的參數(shù)、根據(jù)式(1)和式(2)計(jì)算，據(jù)此可對(duì)n個(gè)樣品發(fā)熱量均值的不確定度進(jìn)行評(píng)定。

(3) 當(dāng)發(fā)熱量樣本不服從正態(tài)分布且n<30時(shí)，要對(duì)n個(gè)樣品發(fā)熱量均值的不確定度進(jìn)行評(píng)定，可采用挑選樣抽法進(jìn)行虛擬采樣，得到n個(gè)樣品發(fā)熱量的均值樣本，然后根據(jù)該均值樣本的概率密度分布對(duì)發(fā)熱量均值的不確定度評(píng)定。

(4) 采用虛擬采樣均值樣本對(duì)n個(gè)樣品發(fā)熱量均值的不確定度評(píng)定結(jié)果與根據(jù)中心極限定理得到的結(jié)論一致。