基于Savitzky-Golay 算法的發(fā)動(dòng)機(jī)動(dòng)態(tài)工況燃燒穩(wěn)定性評價(jià)方法

（寧波吉利羅佑發(fā)動(dòng)機(jī)零部件有限公司 浙江 寧波 315336）

引言

為應(yīng)對全球排放法規(guī)日益加嚴(yán)，各汽車公司對傳統(tǒng)內(nèi)燃機(jī)的瞬時(shí)排放重視程度越來越高，尤其是對發(fā)動(dòng)機(jī)起動(dòng)階段和暖機(jī)階段的排放更為重視。有數(shù)據(jù)表明，WLTC 排放循環(huán)中，50%～80%的一氧化碳[1]排放均來自于發(fā)動(dòng)機(jī)起動(dòng)后包含催化器起燃的30 s內(nèi)。另外，隨著消費(fèi)者對汽車舒適性要求的提高，汽車NVH 也受到越來越多的關(guān)注。發(fā)動(dòng)機(jī)起動(dòng)后的熱機(jī)階段，較大噪聲和振動(dòng)均會(huì)令乘坐者感到不適。以上涉及的問題，均與發(fā)動(dòng)機(jī)燃燒緊密相關(guān)。穩(wěn)定的燃燒，既可以改善排放水平，也對振動(dòng)和噪聲改善明顯。為了評價(jià)燃燒質(zhì)量，本文提出了一種新的發(fā)動(dòng)機(jī)燃燒穩(wěn)定性評價(jià)方法，經(jīng)試驗(yàn)驗(yàn)證，該方法不僅能夠評估發(fā)動(dòng)機(jī)穩(wěn)態(tài)燃燒質(zhì)量，而且還能對動(dòng)態(tài)工況下燃燒穩(wěn)定性做出合理評估，本文定義這種新方法為DIOC，即Dynamic Index of Combustion，下文均以DIOC 代表該評價(jià)方法。

1 目前主要燃燒評價(jià)方法

基于目前燃燒評價(jià)的主要指標(biāo)是IMEPCOV（Coefficient of Variation of Indicated Mean Efficient Pressure，平均指示壓力循環(huán)變動(dòng)率）和IMEPLNV（Least Normalized Value of Indicated Mean Efficient Pressure，平均指示壓力最小歸一化數(shù)值），這兩種評價(jià)指標(biāo)被各汽車主機(jī)廠、研究院和大學(xué)廣泛用于發(fā)動(dòng)機(jī)穩(wěn)態(tài)工況燃燒評價(jià)。IMEPCOV 的計(jì)算公式如式（1）～（3）所示。

其中：σ 代表所測試的N 組平均指示壓力數(shù)據(jù)的標(biāo)準(zhǔn)差，μ 代表所測試的N 組平均指示壓力數(shù)據(jù)的平均值，N 一般取值100～300，yi代表所測N 組數(shù)據(jù)中的任意一組，i=1，2，3，4…N。IMEPLNV 的計(jì)算公式如式（4）所示。

其中：min 代表取最小值運(yùn)算，目的是選取所測試N 組平均指示壓力數(shù)據(jù)中的最小值，μ 代表所測試的N 組平均指示壓力的平均值，N 一般取值100～300。

上述兩種發(fā)動(dòng)機(jī)燃燒評價(jià)方法主要針對發(fā)動(dòng)機(jī)穩(wěn)態(tài)工況進(jìn)行燃燒評價(jià)，無法實(shí)現(xiàn)發(fā)動(dòng)機(jī)動(dòng)態(tài)工況燃燒評價(jià)。這是因?yàn)樯鲜鰞煞N方法均用到了求IMEP平均值的計(jì)算，這意味著評價(jià)時(shí)均需要連續(xù)采集N組試驗(yàn)數(shù)據(jù)，然后對N 組試驗(yàn)數(shù)據(jù)求平均值，得出該穩(wěn)態(tài)工況的期望IMEP 值，根據(jù)實(shí)測IMEP 數(shù)值和期望IMEP 數(shù)值的差值，經(jīng)過計(jì)算后，對燃燒進(jìn)行評價(jià)。對于穩(wěn)態(tài)工況，IMEP 的平均值與期望值數(shù)值相等，但是對于動(dòng)態(tài)工況，IMEP 在一段時(shí)間內(nèi)是連續(xù)大幅變化的數(shù)值，所求平均值會(huì)與真實(shí)IMEP 的期望值不符，這就導(dǎo)致上述兩種方法無法對發(fā)動(dòng)機(jī)動(dòng)態(tài)工況進(jìn)行燃燒評價(jià)。

2 影響發(fā)動(dòng)機(jī)燃燒穩(wěn)定性的因素

發(fā)動(dòng)機(jī)的缸內(nèi)燃燒是一系列復(fù)雜的物理變化和化學(xué)反應(yīng)過程。了解影響發(fā)動(dòng)機(jī)燃燒的關(guān)鍵因素，對于客觀評價(jià)發(fā)動(dòng)機(jī)燃燒和提出優(yōu)化建議十分必要。

2.1 混合氣對發(fā)動(dòng)機(jī)燃燒的影響

影響汽油發(fā)動(dòng)機(jī)燃燒的一個(gè)主要因素是混合氣的空燃比，理論上1 kg 燃料完全燃燒時(shí)需要14.7 kg的空氣。這種空氣和燃料的比例稱為化學(xué)當(dāng)量比。空燃比小于化學(xué)當(dāng)量比時(shí)供給濃混合氣，此時(shí)發(fā)動(dòng)機(jī)發(fā)出的功率大，但燃燒不完全，生成的CO、HC 較多；當(dāng)空燃比大于化學(xué)當(dāng)量比時(shí)供給稀混合氣時(shí)，火焰?zhèn)鞑ニ俣认陆?，燃燒速度變慢，燃燒不穩(wěn)定，混合氣無法完全燃燒，使得HC 增多，甚至發(fā)生失火現(xiàn)象[2]。

2.2 點(diǎn)火性能對發(fā)動(dòng)機(jī)燃燒的影響

點(diǎn)火提前角過小，點(diǎn)火與燃燒在活塞下行且氣缸容積不斷增大的膨脹過程中進(jìn)行。冷機(jī)工況下，混合氣經(jīng)過壓縮到上止點(diǎn)位置時(shí)溫度升高到達(dá)頂峰值，之后溫度會(huì)逐漸降低，導(dǎo)致同樣點(diǎn)火能量下火花塞跳火后火核形成困難，易造成發(fā)動(dòng)機(jī)失火。點(diǎn)火提前角過大，氣缸壓力升高過快，消耗的壓縮功增加，容易誘發(fā)發(fā)動(dòng)機(jī)抖動(dòng)，發(fā)生爆燃的可能性增加，導(dǎo)致發(fā)動(dòng)機(jī)燃燒不穩(wěn)[3]。

2.3 殘余廢氣量對發(fā)動(dòng)機(jī)燃燒的影響

由于燃燒室的存在，發(fā)動(dòng)機(jī)燃燒過程中不可避免地會(huì)有一部分殘余的燃燒廢氣滯留在缸內(nèi)，現(xiàn)代搭載雙VVT 的發(fā)動(dòng)機(jī)，通過更改進(jìn)排氣VVT 重疊角的方法，可以在一定范圍內(nèi)控制滯留在氣缸內(nèi)的殘余廢氣量。冷機(jī)工況下，適量的殘余廢氣可以提升缸內(nèi)溫度，促進(jìn)燃油蒸發(fā)形成均勻混合氣，加快火焰?zhèn)鞑ニ俣龋兄诟纳迫紵€(wěn)定性；過量的殘余廢氣則會(huì)對混合氣形成稀釋效應(yīng)，降低火焰?zhèn)鞑ニ俣?，?dǎo)致燃燒不穩(wěn)定，甚至失火[4]。

3 基于SG 算法的瞬態(tài)工況燃燒穩(wěn)定性評價(jià)

為克服第一節(jié)提出的兩種方法對動(dòng)態(tài)工況燃燒評價(jià)的局限性，特提出基于數(shù)據(jù)擬合方法求IMEP 期望值的方法。數(shù)據(jù)擬合是統(tǒng)計(jì)學(xué)中常用的方法，即根據(jù)實(shí)測IMEP 數(shù)據(jù)，采用基于Savitzky-Golay 算法（以下簡稱SG 算法），對實(shí)測數(shù)據(jù)采用分段多項(xiàng)式擬合，形成一條光滑的IMEP 期望值曲線數(shù)據(jù)。采用基于SG 算法的多項(xiàng)式擬合，可以有效去除設(shè)備測試誤差和傳感器誤差[5]，同時(shí)擬合的IMEP 期望值曲線與實(shí)際IMEP 測試曲線變化趨勢高度吻合。從而實(shí)現(xiàn)動(dòng)態(tài)工況下燃燒質(zhì)量的客觀評價(jià)。這是一種將統(tǒng)計(jì)學(xué)、燃燒和測試3 種學(xué)科結(jié)合的創(chuàng)新思路，瞬態(tài)燃燒穩(wěn)定性評價(jià)的主要步驟如圖1 所示。

圖1 動(dòng)態(tài)工況燃燒評價(jià)流程圖

原始數(shù)據(jù)采集：通過燃燒分析儀采集動(dòng)態(tài)工況原始IMEP 數(shù)值。

SG 擬合：選取合適的參數(shù)，實(shí)現(xiàn)基于SG 算法的最優(yōu)擬合。

擬合質(zhì)量評價(jià)：通過統(tǒng)計(jì)參數(shù)評價(jià)擬合結(jié)果是否滿足工程需要。

動(dòng)態(tài)工況燃燒評價(jià)：動(dòng)態(tài)工況燃燒質(zhì)量是否滿足標(biāo)準(zhǔn)要求。

3.1 原始數(shù)據(jù)采集

本試驗(yàn)研究對象為一臺(tái)1.5L 三缸渦輪增壓氣道噴射發(fā)動(dòng)機(jī)，燃燒數(shù)據(jù)采用奇石樂Kistler5064 D11 系列燃燒分析儀。ECU 數(shù)據(jù)使用INCA 軟件采集，使用MDA 對INCA 采集的數(shù)據(jù)進(jìn)行分析。

3.2 基于Savitzky-Golay 原理的擬合簡介

Savitzky-Golay 算法早年由Savitzky 和Golay 在1964 年發(fā)表，被廣泛地運(yùn)用于數(shù)據(jù)流平滑與除噪，是一種在時(shí)域內(nèi)基于多項(xiàng)式，通過移動(dòng)窗口利用最小二乘法進(jìn)行最佳分段擬合的方法[6]。

3.2.1 SG 擬合理論

通過采集數(shù)據(jù)點(diǎn)yi附近固定個(gè)數(shù)的點(diǎn)擬合一個(gè)多項(xiàng)式，多項(xiàng)式在yi的值就給出了它的光滑數(shù)值μi，即yi點(diǎn)處的期望值（μi∈U 所求光滑的期望值矩陣）。以nl表示yi左邊點(diǎn)的個(gè)數(shù)，以nr表示yi右邊點(diǎn)的個(gè)數(shù)，pi（y）表示相對于點(diǎn)yi的一個(gè)M 次多項(xiàng)式，用最小二乘法擬合nl+nr+1 個(gè)點(diǎn)。因此：

根據(jù)實(shí)測數(shù)據(jù)yi擬合M 次多項(xiàng)式pi（y），可寫成：

用最小二乘法獲得多項(xiàng)式pi（y）來擬合測試數(shù)據(jù)，需要準(zhǔn)確確定bk數(shù)值，使得：

利用這些定義，就可以用矩陣形式進(jìn)行表述此問題：

為求解B，使得式（11）達(dá)到最優(yōu)，令A(yù)B=Y，左乘AT，方程變?yōu)锳TAB=ATY，根據(jù)式（12），則可以進(jìn)一步求得：

經(jīng)過上述計(jì)算，可知矩陣Y=A（ATA）-1ATY，則所求得光滑數(shù)值為:

即為所求光滑數(shù)值矩陣[7]。根據(jù)濾波原理，可得出nl和nr越大，曲線的平滑性越好，不過這是以犧牲曲線特性為代價(jià)的；相反，M 的取值越大，曲線特性失真率越低但使得數(shù)據(jù)低頻部分的平滑度降低，計(jì)算量也變大[8]。目前已經(jīng)有不少成熟軟件集成了SG 算法工具箱，本文采用Uniplot 軟件中的SG 工具箱實(shí)現(xiàn)IMEP 的擬合，如式（14）所示。

其中：μi代表IMEP 期望值，[yi]代表所測N 組數(shù)據(jù)的集合（i=1，2，3，4……N），order 代表多項(xiàng)式擬合的階數(shù)，一般取2 或者3。length 代表分段擬合的窗口寬度，通常取5～10。本文通過對多組實(shí)測數(shù)據(jù)進(jìn)行多個(gè)不同濾波參數(shù)處理后的綜合比較，在滿足確定失真率不大的情況下，盡可能地增加平滑性，最后選擇出了所采用length=5、M=2 的設(shè)定參數(shù)。

3.2.2 SG 擬合與其它動(dòng)態(tài)工況擬合算法比較

目前常用的動(dòng)態(tài)工況擬合算法還有移動(dòng)窗口平均濾波法和一階低通濾波算法，這兩種算法在信號處理和數(shù)據(jù)處理方面應(yīng)用非常廣泛。本文取發(fā)動(dòng)機(jī)動(dòng)態(tài)工況下第2 缸的IMEP 為分析對象，分別采用SG 擬合和另外兩種擬合算法作比較，比較結(jié)果如圖2～圖4 所示。

圖2 移動(dòng)窗口平均濾波結(jié)果

圖3 一階低通濾波結(jié)果

圖4 SG 算法擬合結(jié)果

采用移動(dòng)窗口平均濾波方法的擬合曲線雖然平順性很好，但是在動(dòng)態(tài)工況時(shí)濾波的曲線無法很好跟隨實(shí)際值的變化，存在偏差過大的問題；一階低通濾波相比移動(dòng)窗口濾波方法有所改善，但是同樣存在濾波值滯后實(shí)際測量值的問題。相比較而言，SG擬合的曲線兼具平順性和實(shí)時(shí)性的優(yōu)點(diǎn)，是動(dòng)態(tài)工況燃燒數(shù)據(jù)擬合的最佳方法。

3.3 擬合質(zhì)量評價(jià)方法

采用基于實(shí)測IMEP 值和期望IMEP 值的線性回歸分析法對SG 擬合質(zhì)量進(jìn)行評估。應(yīng)用最小二乘法，線性擬合的通用公式如式（15）所示：

其中：μi代表IMEP 期望值，yi代表所測N 組數(shù)據(jù)的IMEP 實(shí)測值（i=1，2，3，4……N），a 代表回歸直線的斜率，b 代表回歸直線的截距。依據(jù)最優(yōu)相關(guān)系數(shù)r2原則，確定回歸斜率和截距。表1 為線性回歸分析的允差限值。

表1 線性回歸分析的允差限值

3.4 動(dòng)態(tài)工況燃燒評價(jià)指標(biāo)

對于動(dòng)態(tài)工況燃燒質(zhì)量的評價(jià)，本文定義了DIOC 評價(jià)指標(biāo)。動(dòng)態(tài)工況燃燒評價(jià)指標(biāo)DIOC 描述如式（16）所示：

當(dāng)μi≥3 時(shí)，即發(fā)動(dòng)機(jī)平均指示有效壓力期望值≥3 時(shí)，DIOC>90%，認(rèn)定燃燒穩(wěn)定性合格；當(dāng)＜3時(shí)，即發(fā)動(dòng)機(jī)平均指示有效壓力期望值≤3 時(shí)，DIOC>80%，認(rèn)定燃燒穩(wěn)定性合格；當(dāng)發(fā)動(dòng)機(jī)處于催化器起燃工況時(shí)，DIOC>50%認(rèn)定燃燒穩(wěn)定性合格。

4 數(shù)據(jù)分析

圖5 是發(fā)動(dòng)機(jī)在室溫條件下冷態(tài)工況起動(dòng)后的實(shí)際燃燒情況，采用SG 擬合，各項(xiàng)擬合評價(jià)指標(biāo)如表2 所示。

表2 線性回歸分析的實(shí)際值

依據(jù)DIOC 評價(jià)指標(biāo)，圖5 中發(fā)動(dòng)機(jī)在催化器起燃階段存在燃燒穩(wěn)定性低于限值的工況點(diǎn)，需要對這些工況進(jìn)行分析，哪些原因?qū)е氯紵€(wěn)定性變差，然后再制定針對性優(yōu)化方案。

圖5 起動(dòng)暖機(jī)階段動(dòng)態(tài)工況評價(jià)

5 結(jié)束語

本研究針對實(shí)際工作中遇到的發(fā)動(dòng)機(jī)動(dòng)態(tài)工況燃燒穩(wěn)定性評價(jià)問題，采用了基于SG 算法對發(fā)動(dòng)機(jī)實(shí)際IMEP 進(jìn)行擬合，得到IMEP 的期望值曲線，根據(jù)實(shí)測值和期望值的對比，進(jìn)而實(shí)現(xiàn)對動(dòng)態(tài)工況下燃燒穩(wěn)定性的評價(jià)。該方法相比傳統(tǒng)的方法，更能夠客觀評價(jià)動(dòng)態(tài)工況的燃燒穩(wěn)定性，具有較好的效果。雖然本文提出的基于SG 擬合算法的發(fā)動(dòng)機(jī)動(dòng)態(tài)工況燃燒穩(wěn)定性評價(jià)方法尚不完全成熟，需要投入時(shí)間和精力去完善。但是，需要著重強(qiáng)調(diào)的是本文將信號處理領(lǐng)域的思路方法引入到發(fā)動(dòng)機(jī)的燃燒評價(jià)中，這種跨領(lǐng)域的研究嘗試及交叉學(xué)科知識(shí)綜合運(yùn)用，對相關(guān)領(lǐng)域的工作者具有一定的啟發(fā)意義，值得推廣。