朱閃閃，曹永華，謝中敏

（江蘇航空職業(yè)技術(shù)學(xué)院航空工程學(xué)院，江蘇 鎮(zhèn)江 212134）

0 引言

航空發(fā)動(dòng)機(jī)是一種復(fù)雜的裝置，它的研制水平代表了一個(gè)國(guó)家的工業(yè)、科技以及國(guó)防實(shí)力的水平，目前，只有美國(guó)、英國(guó)、中國(guó)等少有的幾個(gè)國(guó)家具有這樣的自主研發(fā)能力[1]。其余各個(gè)國(guó)家的航空界正在努力地尋找一種效率高、成本低、風(fēng)險(xiǎn)低的解決措施。在這一解決措施中，計(jì)算機(jī)仿真具有很大的優(yōu)勢(shì)，通過(guò)計(jì)算機(jī)仿真可以有效地縮短研制周期，減少研制成本[2]。陳玉春等[3]提出了變步長(zhǎng)的牛頓-拉夫遜法對(duì)航空發(fā)動(dòng)機(jī)模型進(jìn)行求解，但是這種方法僅在單軸渦噴發(fā)動(dòng)機(jī)上得到了應(yīng)用，不具備普遍性。施洋等[4]提出罰函數(shù)法等優(yōu)化方法來(lái)代替發(fā)動(dòng)機(jī)非線性方程組的迭代求解過(guò)程，從而完成航空發(fā)動(dòng)機(jī)性能的仿真模擬計(jì)算。這些方法對(duì)計(jì)算初值沒(méi)有要求，在一定程度上緩解了計(jì)算過(guò)程中的發(fā)散問(wèn)題，但是其計(jì)算速度非常緩慢，在實(shí)際應(yīng)用中適用性有限。所以本文首先利用Matlab/Simulink 軟件對(duì)DGEN380 發(fā)動(dòng)機(jī)進(jìn)行建模仿真，然后利用牛頓-拉夫遜方法對(duì)穩(wěn)態(tài)模型進(jìn)行求解，針對(duì)該算法中可能存在結(jié)果不收斂的原因進(jìn)行分析總結(jié)，最后對(duì)穩(wěn)態(tài)求解過(guò)程中的收斂性進(jìn)行仿真分析研究。

1 發(fā)動(dòng)機(jī)模型的建立

目前發(fā)動(dòng)機(jī)的建模大體可以分為機(jī)理建模和實(shí)驗(yàn)建模[5，6]，機(jī)理建模是以發(fā)動(dòng)機(jī)運(yùn)行過(guò)程中所遵循的熱力學(xué)過(guò)程所建立的數(shù)學(xué)模型，而試驗(yàn)建模（系統(tǒng)辨識(shí)）是根據(jù)發(fā)動(dòng)機(jī)的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行分析而獲得的數(shù)學(xué)模型，對(duì)實(shí)驗(yàn)數(shù)據(jù)有著較高的要求。而機(jī)理建模過(guò)程簡(jiǎn)便，所得到的模型精度也較高，并且可以在發(fā)動(dòng)機(jī)運(yùn)行全包線內(nèi)對(duì)發(fā)動(dòng)機(jī)穩(wěn)態(tài)、動(dòng)態(tài)特性進(jìn)行分析。這種方法廣泛地被國(guó)內(nèi)外研究人員所使用。本文將結(jié)合DGEN380 發(fā)動(dòng)機(jī)設(shè)計(jì)點(diǎn)實(shí)驗(yàn)數(shù)據(jù)建立該發(fā)動(dòng)機(jī)的穩(wěn)態(tài)模型。

1.1 發(fā)動(dòng)機(jī)結(jié)構(gòu)與建模假設(shè)

本文的建模對(duì)象是DGEN380 的齒輪傳動(dòng)風(fēng)扇發(fā)動(dòng)機(jī)（GTF），按照其工作原理主要?jiǎng)澐譃檫M(jìn)氣道、風(fēng)扇、齒輪箱[6]、壓氣機(jī)等部件。其結(jié)構(gòu)劃分簡(jiǎn)圖如圖1 所示。

圖1 DGEN380 發(fā)動(dòng)機(jī)結(jié)構(gòu)簡(jiǎn)圖

建立航空發(fā)動(dòng)機(jī)模型時(shí)，作以下假設(shè)：

（1）將發(fā)動(dòng)機(jī)內(nèi)的氣體流動(dòng)視為一維流動(dòng)；

（2）忽略發(fā)動(dòng)機(jī)放氣和引起的影響；

（3）不考慮容積慣性和燃燒延遲的影響；

（4）各個(gè)部件中的比熱值取定值。

1.2 DGEN380 發(fā)動(dòng)機(jī)部件級(jí)模型的建立

（1）初始變量模塊

該模塊主要是由DGEN380 發(fā)動(dòng)機(jī)的高、低壓轉(zhuǎn)速以及四個(gè)旋轉(zhuǎn)部件的壓比組成的初值變量。

（2）減速齒輪箱模塊

根據(jù)理論公式，如果當(dāng)前低壓軸的物理轉(zhuǎn)速為nl，那么風(fēng)扇的物理轉(zhuǎn)速為

（4）流量和效率插值模塊

該模塊包括風(fēng)扇的特性數(shù)據(jù)，因此要根據(jù)當(dāng)前風(fēng)扇的換算轉(zhuǎn)速和壓比來(lái)求得風(fēng)扇當(dāng)前的效率和換算流量[7]，即：

從而根據(jù)式（4）得到風(fēng)扇實(shí)際的效率和流量：

式中k為空氣的絕熱指數(shù)，在此取1.4，R為氣體常數(shù)，在此取287.04 J/kg·K。

在simulink 環(huán)境下建立好的DGEN380 發(fā)動(dòng)機(jī)的部件級(jí)模型如圖2 所示。

圖2 DGEN380 發(fā)動(dòng)機(jī)部件級(jí)模型

2 發(fā)動(dòng)機(jī)穩(wěn)態(tài)模型的建立與求解

2.1 穩(wěn)態(tài)模型的建立

航空發(fā)動(dòng)機(jī)處于穩(wěn)態(tài)時(shí)，需要滿足一定的穩(wěn)態(tài)平衡條件，選取高、低壓轉(zhuǎn)速、風(fēng)扇、壓氣機(jī)、高、低壓渦輪的壓比作為猜值變量，即[nL，πF，nH，πC，πHT，πLT]，故建立以下平衡方程[8]。

式中，ηhpm，ηlpm為高，低壓轉(zhuǎn)子的機(jī)械效率。分別取0.93和0.81。

根據(jù)前面所建立的發(fā)動(dòng)機(jī)整機(jī)模型可知，但飛行條件和燃油流量已知并保持一定時(shí)，只需要知道6 個(gè)狀態(tài)參數(shù)[nL，πF，nH，πC，πHT，πLT]即可得到發(fā)動(dòng)機(jī)當(dāng)前所有的熱力學(xué)參數(shù)，所發(fā)動(dòng)機(jī)的穩(wěn)態(tài)模型可以用這六個(gè)參數(shù)來(lái)表示，即以下非線性方程：

求解上式（10）即是對(duì)發(fā)動(dòng)機(jī)穩(wěn)態(tài)模型求解，從而得到各工作參數(shù)。

2.2 穩(wěn)態(tài)模型的求解

航空發(fā)動(dòng)機(jī)穩(wěn)態(tài)模型求解一般采用數(shù)值法。如N-R法[9]、改進(jìn)的牛頓法[10]等。其中N-R 法最為常用，本文選取N-R 法對(duì)發(fā)動(dòng)機(jī)的穩(wěn)態(tài)模型進(jìn)行求解。

根據(jù)上式（10），令F（x）= [f1，f2，f3，f4，f5，f6]，x= [nL，πF，nH，πC，πHT，πLT]T。N-R 的核心公式為x(k+1)=xk-J-1F（x(k)），其中xk為上一步初猜值的解，為xk+1迭代一步之后修正過(guò)后的解。J為雅可比矩陣，其表達(dá)式為：

采用向前差分法[11]對(duì)各個(gè)變量求偏導(dǎo)，對(duì)初猜解向量中的各項(xiàng)進(jìn)行擾動(dòng)求出雅可比矩陣各列，如第一次擾動(dòng)nL，將擾動(dòng)后的向量[nL+ △nL，πF，nH，πC，πHT，πLT]代入式（11），求出雅可比矩陣的第一列，依次類推，得到最終的雅可比矩陣。

3 穩(wěn)態(tài)求解中可能存在的問(wèn)題分析

3.1 迭代中斷

要保證迭代過(guò)程可以順利地進(jìn)行，首先保證迭代過(guò)程不會(huì)出現(xiàn)中斷，而在仿真的過(guò)程中，經(jīng)常會(huì)出現(xiàn)迭代中斷的問(wèn)題，原因大致可以總結(jié)為以下兩點(diǎn)：

（1）插值過(guò)程無(wú)法進(jìn)行

部件級(jí)模型都是基于部件特性來(lái)建立的，在整機(jī)模型中，風(fēng)扇、壓氣機(jī)、高壓渦輪、低壓渦輪內(nèi)的部件特性的插值能否順利地進(jìn)行，決定著迭代是否完成。

以風(fēng)扇為例，其部件特性的插值如上式（3）所示，如果在插值的過(guò)程中，插值的輸入量超界將會(huì)使插值計(jì)算過(guò)程無(wú)法進(jìn)行，最終導(dǎo)致穩(wěn)態(tài)求解無(wú)法完成。而產(chǎn)生插值中斷的原因主要是由于獨(dú)立變量設(shè)置的不合理，當(dāng)設(shè)置的獨(dú)立變量相比于真實(shí)值偏差較大時(shí)，會(huì)使迭代計(jì)算求解過(guò)程中斷。

（2）迭代過(guò)程中出現(xiàn)奇異

迭代過(guò)程中出現(xiàn)計(jì)算奇異的情況主要是由于開(kāi)負(fù)平方根造成的，例如風(fēng)扇求解流量時(shí)：

如果上式的溫度出現(xiàn)負(fù)值的情況，將會(huì)導(dǎo)致計(jì)算無(wú)法進(jìn)行，迭代中斷。出現(xiàn)這種奇異一般是獨(dú)立變量設(shè)置不當(dāng)造成的。表1 是在發(fā)動(dòng)機(jī)設(shè)計(jì)點(diǎn)處設(shè)置不同的獨(dú)立變量對(duì)計(jì)算結(jié)果的影響。

表1 初值對(duì)計(jì)算結(jié)果的影響

從以上可以看出，合適的獨(dú)立變量是進(jìn)行迭代的第一步，所以對(duì)獨(dú)立變量限制的設(shè)置顯得至關(guān)重要。

3.2 迭代不收斂

大量的實(shí)驗(yàn)數(shù)據(jù)表明，基本的N-R 法雖然完成了迭代的過(guò)程，但是大多數(shù)情況下結(jié)果都是不收斂的，即穩(wěn)態(tài)共同工作方程的求解精度不滿足要求，這樣殘差會(huì)隨著迭代過(guò)程不斷地增大，導(dǎo)致最后的結(jié)果發(fā)散，這一結(jié)果的主要原因有：

（1）獨(dú)立變量初值設(shè)置的不合理

由于發(fā)動(dòng)機(jī)模型的特殊性，在求解迭代過(guò)程中也許會(huì)使求解結(jié)果與現(xiàn)實(shí)狀態(tài)不符，如燃燒室的出口總溫小于進(jìn)口總溫，這時(shí)的結(jié)果只具有數(shù)學(xué)意義而沒(méi)有物理意義的參數(shù)[10]。雖然完成了迭代過(guò)程，但最后結(jié)果是發(fā)散的。

（2）基本N-R 迭代算法本身的不足

N-R 方法的關(guān)鍵是牛頓步長(zhǎng)的選取，一般情況下，牛頓步長(zhǎng)選取是獨(dú)立變量的0.001 倍。而如果獨(dú)立變量的殘差較大時(shí)，擾動(dòng)量較小將起不到什么作用，如果獨(dú)立變量的殘差較小，擾動(dòng)量如果太大將會(huì)導(dǎo)致特性圖插值超界。所以如果采用固定步長(zhǎng)，很容易導(dǎo)致最后結(jié)果的不收斂。

根據(jù)以上所述，導(dǎo)致發(fā)動(dòng)機(jī)穩(wěn)態(tài)求解最后結(jié)果不收斂的情況大致可以歸結(jié)為兩點(diǎn)，變量初值設(shè)置的不合適，二是迭代算法中牛頓步長(zhǎng)選取有問(wèn)題。限制獨(dú)立變量X 的初值的目的是為了保證獨(dú)立變量的初值不會(huì)導(dǎo)致部件插值出現(xiàn)外插值的情況，具體的做法是：例如風(fēng)扇和壓氣機(jī)，根據(jù)特性圖中nccor的范圍，可以對(duì)獨(dú)立變量進(jìn)行限制。對(duì)于高、低壓渦輪，根據(jù)特性圖中Wgcor和ntcor的范圍，可以對(duì)獨(dú)立變量進(jìn)行限制：Wgcor如果超限，則對(duì)其進(jìn)行限制；根據(jù)渦輪進(jìn)口總溫和nc可以計(jì)算得到ntcor，如果ntcor超過(guò)范圍，則限制nc。其他部件的計(jì)算根據(jù)具體情況采取參數(shù)限制措施。獨(dú)立變量初值的限制，保證N-R 法計(jì)算的第一步能夠完成。而對(duì)于牛頓步長(zhǎng)的選取，本文采用變牛頓步長(zhǎng)法使得迭代過(guò)程能夠在保證精度的情況下順利地進(jìn)行。以下根據(jù)這兩個(gè)問(wèn)題進(jìn)行仿真算例分析。

4 仿真算例及結(jié)果分析

根據(jù)以上建立的發(fā)動(dòng)機(jī)穩(wěn)態(tài)模型，對(duì)發(fā)動(dòng)機(jī)的穩(wěn)態(tài)求解過(guò)程進(jìn)行仿真，本文選取飛行包線內(nèi)比較有代表意義的5 個(gè)點(diǎn)，見(jiàn)表2。在基本牛頓-拉夫遜方法的基礎(chǔ)上，采用獨(dú)立變量初值的設(shè)置方法以及牛頓變步長(zhǎng)的方法改善發(fā)動(dòng)機(jī)的穩(wěn)態(tài)求解的迭代過(guò)程。表3 為兩種方法的對(duì)比結(jié)果。

表2 油門桿位置對(duì)應(yīng)的燃油流量

表3 不同方法對(duì)求解結(jié)果的影響

從以上結(jié)果不難看出，對(duì)獨(dú)立變量的限制可以加快迭代速度，改善迭代效果。迭代算法的改進(jìn)（即在迭代過(guò)程中采用變步長(zhǎng)），可以使迭代朝著更加有效的方向進(jìn)行。所以以上方法不僅有效地改善了收斂性，還對(duì)收斂速度有了進(jìn)一步的改善。

5 結(jié)論

本文對(duì)DGEN380 發(fā)動(dòng)機(jī)部件特性計(jì)算中的不收斂問(wèn)題進(jìn)行了系統(tǒng)的研究后，以牛頓- 拉夫遜方法為基礎(chǔ)，對(duì)該方法進(jìn)行了改進(jìn)并取得了非常好的效果，改善了算法迭代收斂性的問(wèn)題。認(rèn)為可以不必使用其它計(jì)算量大的求解算法。同樣這些改進(jìn)辦法也可以為其他類型的航空發(fā)動(dòng)機(jī)的特性計(jì)算提供參考。