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

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

引言

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

1 目前主要燃燒評價方法

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

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

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

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

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

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

2.1 混合氣對發(fā)動機燃燒的影響

影響汽油發(fā)動機燃燒的一個主要因素是混合氣的空燃比，理論上1 kg 燃料完全燃燒時需要14.7 kg的空氣。這種空氣和燃料的比例稱為化學當量比。空燃比小于化學當量比時供給濃混合氣，此時發(fā)動機發(fā)出的功率大，但燃燒不完全，生成的CO、HC 較多；當空燃比大于化學當量比時供給稀混合氣時，火焰?zhèn)鞑ニ俣认陆担紵俣茸兟紵环€(wěn)定，混合氣無法完全燃燒，使得HC 增多，甚至發(fā)生失火現(xiàn)象[2]。

2.2 點火性能對發(fā)動機燃燒的影響

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

2.3 殘余廢氣量對發(fā)動機燃燒的影響

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

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

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

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

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

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

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

動態(tài)工況燃燒評價：動態(tài)工況燃燒質(zhì)量是否滿足標準要求。

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

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

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

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

3.2.1 SG 擬合理論

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

3.4 動態(tài)工況燃燒評價指標

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

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

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

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

表2 線性回歸分析的實際值

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

圖5 起動暖機階段動態(tài)工況評價

5 結(jié)束語

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