曹韻奇，劉桂娟，段天英

(中國(guó)原子能科學(xué)研究院 反應(yīng)堆工程技術(shù)研究部，北京 102413)

由于鈉冷快堆設(shè)計(jì)的獨(dú)特性，當(dāng)蒸汽發(fā)生器內(nèi)部傳熱管路出現(xiàn)破損時(shí)，將發(fā)生劇烈的鈉-水反應(yīng)[1-2]。其振動(dòng)信號(hào)中除鈉-水反應(yīng)信號(hào)外，還伴隨著具有較高能量的背景噪聲信號(hào)。已有研究[3-4]表明，該背景噪聲信號(hào)由多種聲源交混構(gòu)成，頻譜復(fù)雜。工作溫度、壓力以及流速等的變化均會(huì)對(duì)背景噪聲信號(hào)特征產(chǎn)生極大影響，因此利用傳統(tǒng)的時(shí)域或頻域分析方法很難在復(fù)雜背景噪聲存在條件下快速識(shí)別鈉-水反應(yīng)信號(hào)，從而導(dǎo)致泄漏探測(cè)出現(xiàn)延遲，造成泄漏擴(kuò)散，最終導(dǎo)致事故蔓延。因此，如果能通過(guò)對(duì)蒸汽發(fā)生器運(yùn)行過(guò)程中所產(chǎn)生的信號(hào)進(jìn)行分析并對(duì)后續(xù)信號(hào)可能具有的信號(hào)特征進(jìn)行預(yù)測(cè)，對(duì)于正確分析鈉-水反應(yīng)能量傳播機(jī)理和進(jìn)行蒸汽發(fā)生器快速泄漏判斷具有非常重要的意義。時(shí)間序列建模預(yù)測(cè)分析主要包括分類(lèi)、聚類(lèi)和回歸預(yù)測(cè)等，是基于隨機(jī)過(guò)程理論和數(shù)理統(tǒng)計(jì)學(xué)的一種信號(hào)分析方法，該方法通過(guò)對(duì)觀測(cè)得到的數(shù)據(jù)序列進(jìn)行建模，從動(dòng)態(tài)的角度去發(fā)現(xiàn)序列的變化規(guī)律及趨勢(shì)，從而對(duì)觀測(cè)序列的未來(lái)行為進(jìn)行預(yù)測(cè)。目前時(shí)間序列預(yù)測(cè)分析在氣象預(yù)測(cè)、工程故障診斷、金融股票走勢(shì)分析等數(shù)據(jù)隨時(shí)間變化快且具有一定統(tǒng)計(jì)學(xué)規(guī)律的快變趨勢(shì)預(yù)測(cè)方面具有非常廣泛的應(yīng)用[5-7]。

為解決目前蒸汽發(fā)生器鈉-水反應(yīng)噪聲探測(cè)方法存在的觀測(cè)信號(hào)受背景噪聲干擾嚴(yán)重從而導(dǎo)致泄漏探測(cè)響應(yīng)速度慢、泄漏信號(hào)具體響應(yīng)時(shí)間無(wú)法準(zhǔn)確定位等問(wèn)題，本文將基于探測(cè)信號(hào)的時(shí)間序列建模預(yù)測(cè)分析方法引入鈉-水反應(yīng)噪聲探測(cè)信號(hào)分析。根據(jù)已知觀測(cè)信號(hào)序列進(jìn)行建模來(lái)預(yù)測(cè)未來(lái)信號(hào)的發(fā)展趨勢(shì)，從而對(duì)蒸汽發(fā)生器內(nèi)部是否發(fā)生泄漏進(jìn)行判斷，提高噪聲探測(cè)信號(hào)分析的響應(yīng)速度，識(shí)別、定位泄漏信號(hào)的具體響應(yīng)時(shí)間。

1 時(shí)間序列預(yù)測(cè)分析基本原理

1.1 時(shí)間序列分析的基本形式

時(shí)間序列是指隨時(shí)間改變而發(fā)生隨機(jī)變化的序列，對(duì)于時(shí)間序列分析首要目的是針對(duì)觀測(cè)數(shù)據(jù)序列建立最能體現(xiàn)時(shí)間序列規(guī)律的數(shù)學(xué)模型[8-9]。

對(duì)于任何均值為零的平穩(wěn)時(shí)間序列{Wt}，若序列{Wt}可表示為過(guò)去p時(shí)刻上數(shù)值Wt-1，Wt-2，…，Wt-p以及t時(shí)刻及t之前q+1個(gè)時(shí)刻白噪聲at，at-1，…，at-p的線性組合，則有：

Wt-φ1Wt-1-φ2Wt-2-…-φpWt-p=
at-θ1at-1-θ2at-2-…-θqat-q
t=0，±1，±2，…

(1)

其中：φi(i=1，2，…，n)為模型的自回歸系數(shù)；θi(i=1，2，…，n)為模型的滑動(dòng)平均系數(shù)。

式(1)稱(chēng)為自回歸滑動(dòng)平均模型(ARMA模型)，特別是當(dāng)模型的滑動(dòng)平均系數(shù)項(xiàng)為0時(shí)稱(chēng)為自回歸模型(AR模型)，當(dāng)模型的自回歸系數(shù)項(xiàng)為0時(shí)稱(chēng)為滑動(dòng)平均模型(MA模型)。ARMA模型是目前最為常用的一種時(shí)間序列預(yù)測(cè)模型，該模型適用于隨機(jī)特征不隨時(shí)間變化的平穩(wěn)、單自變量時(shí)間序列建模[10]，但對(duì)于非平穩(wěn)時(shí)間序列的預(yù)測(cè)問(wèn)題往往會(huì)出現(xiàn)統(tǒng)計(jì)特征離散導(dǎo)致的較大預(yù)測(cè)誤差。對(duì)鈉-水反應(yīng)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行單位根檢驗(yàn)(unit root test)后發(fā)現(xiàn)實(shí)驗(yàn)數(shù)據(jù)符合平穩(wěn)時(shí)間序列假設(shè)，因此本研究采用ARMA建模預(yù)測(cè)方法。

1.2 模型的識(shí)別與定階

對(duì)時(shí)間序列進(jìn)行建模分析的第1步是進(jìn)行模型識(shí)別，即在3種常用模型中選取合適的模型類(lèi)型。模型識(shí)別的主要依據(jù)是時(shí)間序列的樣本自相關(guān)函數(shù)(ACF)和樣本偏相關(guān)函數(shù)(PACF)。判斷各模型適用性的ACF和PACF統(tǒng)計(jì)特征列于表1。

表1 ARMA模型的ACF和PACF統(tǒng)計(jì)特征Table 1 Statistical characteristics of ACF and PACF in ARMA model

確定適合的時(shí)間序列模型類(lèi)型后，還要計(jì)算適用的模型階數(shù)，即確定p、q的大小，實(shí)際上，當(dāng)數(shù)據(jù)長(zhǎng)度較大時(shí)根據(jù)自/偏相關(guān)函數(shù)來(lái)確定模型階數(shù)效果并不理想，這是因?yàn)楫?dāng)數(shù)據(jù)長(zhǎng)度較大時(shí)相關(guān)函數(shù)將會(huì)在后端出現(xiàn)振蕩，影響對(duì)拖尾或截尾性質(zhì)的判斷。目前模型階次的常用判斷方法主要有最終預(yù)測(cè)誤差(FPE)準(zhǔn)則、AIC(Akaike information criterion)準(zhǔn)則等方法[11-12]，F(xiàn)PE準(zhǔn)則函數(shù)和AIC準(zhǔn)則函數(shù)表達(dá)式分別為：

(2)

CAIC(m)=Nlnσ2+2m

(3)

式中：N為數(shù)據(jù)長(zhǎng)度；σ2為N的方差。根據(jù)信號(hào)檢驗(yàn)原理，當(dāng)CFPE(m)或CAIC(m)取最小值時(shí)所對(duì)應(yīng)的判斷數(shù)據(jù)長(zhǎng)度m為適用階次[13]。為保證結(jié)果準(zhǔn)確，本文將同時(shí)使用AIC和FPE兩種準(zhǔn)則對(duì)模型進(jìn)行階次判斷。

1.3 時(shí)間序列預(yù)測(cè)

Wt=at+φ1Wt-1+φ2Wt-2+…+φpWt-p-

θ1at-1-θ2at-2-…-θqat-q

(4)

令t=k+s，將第k+s時(shí)刻的白噪聲估計(jì)值ak+s=0代入式(4)，則有：

θ1ak+s-1-θ2ak+s-2-…-θqak+s-q

(5)

θ1ak-θ2ak-1-…-θqak+1-q

(6)

2 實(shí)驗(yàn)分析

2.1 建模預(yù)測(cè)思路

首先對(duì)實(shí)驗(yàn)獲得的蒸汽發(fā)生器原始觀測(cè)信號(hào)進(jìn)行截取后計(jì)算相關(guān)系數(shù)，確定提取觀測(cè)信號(hào)適用的模型及參數(shù)，這樣就可根據(jù)已有的觀測(cè)序列得到適用模型，然后通過(guò)時(shí)間序列預(yù)測(cè)方法計(jì)算后續(xù)時(shí)間信號(hào)的預(yù)測(cè)值并與后續(xù)信號(hào)的實(shí)際觀測(cè)值進(jìn)行比較。實(shí)際上由于實(shí)驗(yàn)開(kāi)始時(shí)并沒(méi)有進(jìn)行注水模擬泄漏，因此觀測(cè)信號(hào)的最初段為該工況下的純背景噪聲信號(hào)，通過(guò)序列建模預(yù)測(cè)得到的預(yù)測(cè)值與下一時(shí)刻觀測(cè)值進(jìn)行對(duì)比，如觀測(cè)值在預(yù)測(cè)范圍之內(nèi)則可認(rèn)為該時(shí)刻未發(fā)生泄漏，此時(shí)刻的觀測(cè)信號(hào)仍是背景噪聲信號(hào)；如觀測(cè)值不在預(yù)測(cè)范圍之內(nèi)則該時(shí)刻引入了新的振動(dòng)信號(hào)源，即泄漏發(fā)生，該時(shí)刻為探測(cè)到泄漏發(fā)生的起始時(shí)刻。對(duì)實(shí)驗(yàn)信號(hào)進(jìn)行分析，結(jié)果表明,應(yīng)用該方法能極大程度提高強(qiáng)背景噪聲存在條件下的鈉-水反應(yīng)泄漏信號(hào)識(shí)別速度，在實(shí)驗(yàn)條件下探測(cè)響應(yīng)時(shí)間小于1 s。該方法步驟如下。

1) 提取未進(jìn)行泄漏模擬的原始觀測(cè)振動(dòng)信號(hào)為樣本序列(為保證實(shí)驗(yàn)數(shù)據(jù)穩(wěn)定，選取實(shí)驗(yàn)開(kāi)始后1～2 s內(nèi)數(shù)據(jù)作為樣本序列)，并對(duì)觀測(cè)樣本信號(hào)采用現(xiàn)代模型法進(jìn)行建模。

2) 根據(jù)樣本序列段建立的模型對(duì)下一時(shí)刻的待預(yù)測(cè)信號(hào)進(jìn)行模型預(yù)測(cè)，得到預(yù)測(cè)信號(hào)。

3) 將預(yù)測(cè)信號(hào)與下一時(shí)刻的觀測(cè)信號(hào)進(jìn)行對(duì)比，如觀測(cè)振動(dòng)信號(hào)在預(yù)測(cè)允許范圍內(nèi)，則將待預(yù)測(cè)的點(diǎn)在觀測(cè)振動(dòng)信號(hào)中向后移動(dòng)一位得到新的待預(yù)測(cè)信號(hào)，進(jìn)行步驟4；如超出預(yù)測(cè)允許范圍則停止計(jì)算。

4) 重復(fù)步驟2、3，直到出現(xiàn)觀測(cè)振動(dòng)信號(hào)不在預(yù)測(cè)范圍內(nèi)則停止計(jì)算并輸出該數(shù)據(jù)點(diǎn)所在的觀測(cè)序列位置，該時(shí)刻為識(shí)別到的泄漏起始時(shí)刻。

2.2 數(shù)據(jù)的獲取及預(yù)處理

為最大程度上模擬快堆蒸汽發(fā)生器運(yùn)行環(huán)境，采用清華大學(xué)鈉-水反應(yīng)實(shí)驗(yàn)臺(tái)架模擬最大背景噪聲存在條件下的鈉-水反應(yīng)觀測(cè)數(shù)據(jù)進(jìn)行實(shí)驗(yàn)。實(shí)驗(yàn)開(kāi)始時(shí)風(fēng)機(jī)等機(jī)械設(shè)施全部開(kāi)啟，鈉流速為6.1 m3/h，水流速為72 L/h，回路鈉溫為250 ℃，此工況能為實(shí)驗(yàn)臺(tái)架提供最大的運(yùn)行背景噪聲，噪聲實(shí)驗(yàn)開(kāi)始后約2.5 s通過(guò)直徑為3 mm的毛細(xì)管向鈉中注水模擬泄漏發(fā)生，注水流量為400 mL/h。實(shí)驗(yàn)系統(tǒng)采樣頻率為3 kHz、采樣時(shí)間為10 s、采樣長(zhǎng)度為30 000數(shù)據(jù)點(diǎn)，原始實(shí)驗(yàn)數(shù)據(jù)經(jīng)零均值處理后如圖1所示。

圖1 最大背景噪聲條件下的原始泄漏信號(hào)Fig.1 Original leakage signal under maximum background noise

2.3 模型的建立

首先確定觀測(cè)序列的適用模型，提取觀測(cè)序列數(shù)據(jù)3 000點(diǎn)作為初始樣本。已有研究[14-15]表明，當(dāng)樣本長(zhǎng)度較大時(shí)，模型階次通常不會(huì)超過(guò)觀測(cè)序列長(zhǎng)度的1/10，對(duì)初始序列進(jìn)行300階次自相關(guān)和偏相關(guān)函數(shù)計(jì)算，結(jié)果如圖2所示。

圖2 樣本序列自相關(guān)和偏相關(guān)函數(shù)Fig.2 Auto-correlation and partial-correlation functions of sample sequence

根據(jù)圖2可確定，樣本序列存在明顯的自相關(guān)函數(shù)拖尾、偏相關(guān)函數(shù)截尾趨勢(shì)，因此對(duì)觀測(cè)序列采用AR模型進(jìn)行建模，分別根據(jù)FPE準(zhǔn)則和AIC準(zhǔn)則進(jìn)行適用階次分析，結(jié)果如圖3所示。

圖3 FPE和AIC準(zhǔn)則函數(shù)變化曲線Fig.3 FPE and AIC criterion function curves

兩種計(jì)算準(zhǔn)則給出的最小值均出現(xiàn)在m=33處，因此觀測(cè)序列的模型階次確定為33。

2.4 模型的預(yù)測(cè)

建模完成后，通過(guò)計(jì)算可得到模型參數(shù)，根據(jù)前文實(shí)驗(yàn)步驟采用單步時(shí)間序列預(yù)測(cè)法進(jìn)行序列預(yù)測(cè)。通過(guò)計(jì)算發(fā)現(xiàn)，當(dāng)數(shù)據(jù)點(diǎn)為6 954時(shí)預(yù)測(cè)值與實(shí)際觀測(cè)值差值不在模型白噪聲范圍內(nèi)，即系統(tǒng)在第6 954數(shù)據(jù)點(diǎn)識(shí)別到有新的聲源加入，蒸汽發(fā)生器內(nèi)部出現(xiàn)鈉-水反應(yīng)。由采樣率Fs=3 000可知，預(yù)測(cè)序列識(shí)別的泄漏時(shí)刻為6 954/3 000≈2.31 s，即實(shí)驗(yàn)開(kāi)始后2.31 s識(shí)別到鈉-水反應(yīng)，響應(yīng)時(shí)間約為0.2 s。模型預(yù)測(cè)法在計(jì)算停止前得到的預(yù)測(cè)信號(hào)如圖4所示。

圖4 基于時(shí)間序列的預(yù)測(cè)信號(hào)Fig.4 Forecast signal based on time series

3 對(duì)比及討論

3.1 處理結(jié)果分析

圖5為實(shí)驗(yàn)開(kāi)始后2～3 s的預(yù)測(cè)序列與原始觀測(cè)序列。通過(guò)計(jì)算得到預(yù)測(cè)誤差的最大值為9.09%，由此可認(rèn)為利用所建立的模型可較好地預(yù)測(cè)信號(hào)趨勢(shì)，預(yù)測(cè)信號(hào)具有良好的置信度，所建立的模型適用于信號(hào)的預(yù)測(cè)分析。

圖5 預(yù)測(cè)序列與原始觀測(cè)序列Fig.5 Forecast sequence and original observation sequence

觀察圖5可發(fā)現(xiàn)，觀測(cè)序列相比于預(yù)測(cè)序列更加離散，預(yù)測(cè)序列能明顯體現(xiàn)出信號(hào)特征變化趨勢(shì)。通過(guò)計(jì)算得到原始觀測(cè)序列的標(biāo)準(zhǔn)差為49.112 7、方差為241.2，預(yù)測(cè)序列的標(biāo)準(zhǔn)差為45.193 6、方差為204.2，說(shuō)明預(yù)測(cè)序列相比于觀測(cè)序列更加穩(wěn)定，接近于統(tǒng)計(jì)期望，因此預(yù)測(cè)序列對(duì)于新的聲源引入更加敏感，在泄漏發(fā)生時(shí)系統(tǒng)能更快給出響應(yīng)。

觀測(cè)序列2～3 s內(nèi)的均方根(RMS)曲線如圖6所示。

圖6 觀測(cè)信號(hào)RMS曲線Fig.6 Observation signal RMS curve

加速度傳感器所測(cè)信號(hào)的RMS曲線與聲壓信號(hào)能量表示形式相同，因此可通過(guò)RMS曲線近似表示信號(hào)的能量強(qiáng)度變化。從圖6可發(fā)現(xiàn)，觀測(cè)信號(hào)能量在2.3 s附近明顯增強(qiáng)，這與采用預(yù)測(cè)序列得到的結(jié)果相同，證明了噪聲實(shí)驗(yàn)基于時(shí)間序列預(yù)測(cè)分析的準(zhǔn)確性和可靠性。并且由于該方法只需對(duì)已有測(cè)量數(shù)據(jù)進(jìn)行建模，即可對(duì)未來(lái)信號(hào)趨勢(shì)進(jìn)行預(yù)測(cè)。綜上可見(jiàn)，該方法在鈉-水反應(yīng)故障信號(hào)的在線檢測(cè)方面具有一定優(yōu)勢(shì)。

3.2 處理結(jié)果對(duì)比

目前對(duì)于鈉-水反應(yīng)噪聲探測(cè)信號(hào)數(shù)據(jù)處理普遍采用小波分解重構(gòu)方法，小波分解方法具有能有效分離不同頻率信號(hào)且同時(shí)具有時(shí)頻域分辨率的特點(diǎn)，因此在工程故障診斷、信號(hào)消噪等方面廣泛應(yīng)用，是目前使用最為廣泛的信號(hào)處理方法[16]。

對(duì)觀測(cè)數(shù)據(jù)通過(guò)Daubechies(DB)小波進(jìn)行一次分解重構(gòu)，得到2～3 s內(nèi)原始觀測(cè)信號(hào)的高頻部分重構(gòu)信號(hào)，如圖7所示。

通過(guò)計(jì)算得到基于Daubechies小波的高頻分量重構(gòu)信號(hào)標(biāo)準(zhǔn)差為17.254 4、方差為297.714 9，與預(yù)測(cè)信號(hào)序列相比，小波分解重構(gòu)信號(hào)各數(shù)據(jù)點(diǎn)更集中，但數(shù)據(jù)整體與統(tǒng)計(jì)均值的偏離更大，因此整體時(shí)域信號(hào)表現(xiàn)更加發(fā)散。預(yù)測(cè)信號(hào)序列和小波重構(gòu)信號(hào)序列的RMS譜如圖8所示。從圖8可發(fā)現(xiàn)，觀測(cè)信號(hào)經(jīng)小波重構(gòu)后能量變化雖更加明顯，但離散程度大，無(wú)法判斷信號(hào)能量突變的具體時(shí)間；預(yù)測(cè)信號(hào)能量變化相較小波重構(gòu)信號(hào)更加穩(wěn)定，同時(shí)在2.3 s附近有明顯的能量突變，說(shuō)明時(shí)間序列預(yù)測(cè)方法得到的信號(hào)能量更加穩(wěn)定，具有良好的預(yù)測(cè)精度和較好的去噪效果。

圖7 高頻分量重構(gòu)信號(hào)Fig.7 High frequency component reconstruction signal

圖8 預(yù)測(cè)信號(hào)與小波重構(gòu)信號(hào)的RMS譜Fig.8 RMS spectra of predictive signal and wavelet reconstruction signal

4 結(jié)論

本文使用時(shí)間序列建模方法對(duì)鈉-水反應(yīng)信號(hào)進(jìn)行了預(yù)測(cè)分析，系統(tǒng)能在強(qiáng)背景噪聲存在條件下對(duì)信號(hào)能量變化進(jìn)行快速響應(yīng)，預(yù)測(cè)序列與觀測(cè)信號(hào)序列相比具有更好的穩(wěn)定性，能更準(zhǔn)確地反映信號(hào)能量的變化趨勢(shì)，對(duì)新的聲源引入更加敏感。相較于傳統(tǒng)時(shí)域判斷方法，應(yīng)用該方法對(duì)信號(hào)進(jìn)行預(yù)測(cè)分析能獲得更好的響應(yīng)速度，并實(shí)現(xiàn)鈉-水反應(yīng)信號(hào)在線檢測(cè)，快速識(shí)別并準(zhǔn)確判斷蒸汽發(fā)生器內(nèi)部泄漏，實(shí)驗(yàn)響應(yīng)時(shí)間小于1 s，為快堆蒸汽發(fā)生器泄漏探測(cè)提供了一種快速識(shí)別方法，也是時(shí)間序列預(yù)測(cè)方法應(yīng)用方面的一種全新嘗試。