趙冬冬，趙國勝，夏磊，方淳，馬睿，皇甫宜耿

1.西北工業(yè)大學(xué) 自動化學(xué)院，西安 710129 2.航空工業(yè)第一飛機設(shè)計研究院，西安 710089

燃料電池是一種通過化學(xué)反應(yīng)，將化學(xué)能轉(zhuǎn)變?yōu)殡娔艿男滦桶l(fā)電裝置[1-2]，目前研究的燃料電池有固體氧化物燃料電池、熔融碳酸鹽燃料電池和質(zhì)子交換膜燃料電池(PEMFC)等多個類型[3]，其中，PEMFC因其效率高、噪聲小、排放無污染等優(yōu)點，被認(rèn)為是未來最具有潛力的無人機用動力源，成為當(dāng)前的研究熱點[4-6]。PEMFC系統(tǒng)通常包含空氣供應(yīng)子系統(tǒng)、氫氣供應(yīng)子系統(tǒng)、濕度控制子系統(tǒng)、溫度控制子系統(tǒng)和功率調(diào)節(jié)子系統(tǒng)5部分。其中，陰極空氣供應(yīng)子系統(tǒng)由空壓機、驅(qū)動電機、供氣歧管、返回歧管等多個部分構(gòu)成[7]，空壓機提供的氧氣流量及陰極氧氣分壓是決定PEMFC系統(tǒng)效率和動態(tài)性能的關(guān)鍵因素[8]。

無人機飛行過程中，由高度變化引起的氣壓變化、溫度變化、濕度變化和氧含量變化等，會直接影響到空壓機的工作特性，進而對燃料電池的工作性能產(chǎn)生影響[9-10]。另一方面，無人機飛行狀態(tài)的變化會直接導(dǎo)致燃料電池負(fù)載的變化，在燃料電池的加載過程中，能否及時根據(jù)負(fù)載的變化對陰極供給流量以及陰極壓力做出調(diào)整，將會直接影響到燃料電池的工作性能[4]。因此，要保證燃料電池系統(tǒng)的高效率和可靠性，對變高度下空壓機的特性以及陰極供氣系統(tǒng)的研究具有重要意義。

空壓機模型一般可通過對大量不同工況下的空壓機實驗數(shù)據(jù)擬合得到，但由于空壓機試驗條件的限制，在全運行范圍內(nèi)的實驗數(shù)據(jù)往往非常有限。此外，查找表的插值程序不是連續(xù)可微的，尤其在實驗數(shù)據(jù)有限的情況下，外推不可靠，查找表的預(yù)測性會大大降低。Deng等[11]基于數(shù)據(jù)驅(qū)動建立了外源輸入的非線性自回歸移動平均(NARMAX)空壓機模型，并采用循環(huán)神經(jīng)網(wǎng)絡(luò)(RNN)方法估計模型參數(shù)。但空壓機進口工況存在很強的相關(guān)性，入口溫度與壓力的變化會影響入口氣體密度和流量的大小，增加了在全轉(zhuǎn)速、變高度下神經(jīng)元的訓(xùn)練難度。CFX、Fluent以及Numeca等軟件通過計算流體動力學(xué)(CFD)的數(shù)值模擬方法模擬空壓機輸出性能，但該方法計算量大，結(jié)果可信度和計算精度都有待提高[12]。偏最小二乘法(PLS)回歸算法是一種由化學(xué)領(lǐng)域發(fā)展形成的新型數(shù)據(jù)統(tǒng)計分析方法，其在預(yù)測空壓機性能方面具有強魯棒性，需要的數(shù)據(jù)量少且具有明顯物理意義的優(yōu)點。Li等[13]在傳統(tǒng)PLS模型的基礎(chǔ)上，使用冪函數(shù)多項式或三角函數(shù)多項式作為基函數(shù)(PLSN)，預(yù)測離心空壓機的壓比和效率，PLSN在計算時間最短的情況下預(yù)測精度最高，模型結(jié)構(gòu)也更加簡單。但PLSN方法并不是對所有的多重共線性問題都有效，某些數(shù)據(jù)處理的結(jié)果會對空壓機的建模影響較大。Gravdahl[14]和Chu[15]等在地面常溫常壓條件下，基于熱力學(xué)第一定律，從能量傳遞的角度提出了一種能夠在虛擬實驗臺(VTB)計算環(huán)境中進行系統(tǒng)仿真的離心式空壓機動態(tài)模型。本文基于文獻[14-15]中的方法預(yù)測摩擦損失，同時考慮了不同高度下溫度、氣壓、密度、雷諾數(shù)等環(huán)境因素變化，建立了適用于不同高度下的空壓機模型。

無人機飛行狀態(tài)的改變會導(dǎo)致燃料電池負(fù)載的劇烈變化，當(dāng)負(fù)載突變時，如果不能根據(jù)負(fù)載電流的變化，及時對進入陰極的空氣流量做出調(diào)整，則會出現(xiàn)氧饑餓或氧飽和現(xiàn)象。氧饑餓現(xiàn)象的發(fā)生會導(dǎo)致燃料電池輸出電壓及功率的降低，嚴(yán)重時更會影響燃料電池的壽命[16]。氧飽和現(xiàn)象的發(fā)生則會導(dǎo)致功率的浪費，進而降低燃料電池系統(tǒng)的效率[17]。為避免氧饑餓和氧飽和現(xiàn)象的發(fā)生，一些文章中提出了根據(jù)負(fù)載變化對PEMFC陰極過氧比(空壓機提供的氧氣量與燃料電池實際反應(yīng)所需的氧氣量之比)的控制。Zhang[18]和Baroud[19]等給出了PEMFC陰極供氣系統(tǒng)的四階模型，并基于所給模型，針對過氧比的控制提出了不同的模糊PID控制算法。劉秋秀[20]采用半機理半實驗建模的方法給出了PEMFC九階供氣系統(tǒng)模型，并提出了基于系統(tǒng)傳遞函數(shù)的一階最優(yōu)自抗擾控制(ADRC)控制策略。王帥[21]給出了四階供氣系統(tǒng)模型，并考慮了空壓機電壓不能無限制增加的因素，提出了具有約束處理能力的模型預(yù)測控制(MPC)。李克雷[22]和郭愛[23]等為保持燃料電池系統(tǒng)凈輸出功率最大，提出了針對過氧比的電堆電流前饋控制和空氣流量模糊反饋控制的策略。浙江大學(xué)Chen等[3]給出了陰極供氣系統(tǒng)的三階模型，并基于反饋線性化的方法設(shè)計了一個非線性控制器來跟蹤控制過氧比保持在最優(yōu)值。Zhang等[24]針對燃料電池供氣系統(tǒng)過氧比提出了一種基于2型模糊邏輯系統(tǒng)(T2-FLS)的自適應(yīng)魯棒控制器，基于李雅普諾夫理論分析了該控制系統(tǒng)的穩(wěn)定性，并通過實驗驗證了該控制器的實用性和可行性。方思雨[25]提出了一種基于狀態(tài)量估計器的準(zhǔn)無限時域模型預(yù)測控制(QHMPC)，采用狀態(tài)量估計器有效解決了QHMPC算法中一些狀態(tài)量不易觀測的問題，并針對節(jié)氣門控制設(shè)計了相應(yīng)控制器。Gruber等[26]針對過氧比的控制提出了一種線性模型預(yù)測控制算法，所提出的線性MPC采用參數(shù)自適應(yīng)的方法來適應(yīng)電堆電流的變化。綜上所述，目前針對PEMFC陰極供氣系統(tǒng)的研究主要集中在對過氧比的控制，對陰極壓力控制的研究較少。然而在無人機飛行過程中，高度的變化會直接引起陰極氣壓的波動，進而對燃料電池輸出性能及功率產(chǎn)生較大影響。因此，除了對過氧比的控制，陰極壓力的控制也是保證無人機用PEMFC系統(tǒng)高效運行的關(guān)鍵因素[27-28]。

本文首先根據(jù)離心空壓機內(nèi)氣體流動特性，通過預(yù)測的總熵增和損失建立了不同高度下的空壓機MAP圖。其次，基于無刷直流電機反電勢特征構(gòu)建了高速空壓機驅(qū)動電機模型，并將空壓機模型與驅(qū)動電機模型結(jié)合至陰極空氣供應(yīng)系統(tǒng)模型中。最后，所提出的陰極供氣系統(tǒng)FFM控制可實現(xiàn)PEMFC在0～3 km寬工況條件下運行，滿足不同飛行狀態(tài)下的功率需求，保證在無人機飛行狀態(tài)和高度發(fā)生改變時，PEMFC過氧比和陰極壓力能夠在短時間內(nèi)跟蹤上給定值，減小由氣壓變化以及供氣不足等對燃料電池壽命和耐久性的影響。

1 無人機用PEMFC陰極供氣系統(tǒng)建模

針對無人機飛行高度的變化，建立了高速離心空壓機的數(shù)學(xué)模型，并給出了不同高度下的空壓機特性分析，根據(jù)能量守恒和質(zhì)量連續(xù)方程，對供應(yīng)歧管、燃料電池電堆和返回歧管進行了動態(tài)建模，具體PEMFC供氣系統(tǒng)框圖如圖1所示。圖中：Wcp為空壓機流量；Wsm.out為供應(yīng)歧管出口流量；Wca.out為陰極出口流量；Wrm.out為返回歧管出口流量;Wan.in和Wan.out分別為陽極入口和出口流量。

圖1 PEMFC系統(tǒng)框圖

1.1 變高度下空壓機與驅(qū)動電機建模

目前針對PEMFC陰極供氣系統(tǒng)中采用的空壓機模型，主要是通過神經(jīng)網(wǎng)絡(luò)等數(shù)據(jù)擬合方法得到，其建模過程存在數(shù)據(jù)量大、處理難度高等缺點，建立的空壓機模型大多只適用于地面，未考慮環(huán)境變化對其帶來的影響。因此，建立變高度下空壓機及驅(qū)動電機模型并進行控制，對燃料電池系統(tǒng)在中低空環(huán)境下的實際應(yīng)用具有重要意義。

1.1.1 環(huán)境模型

空壓機入口氣體特性隨著高度的變化而改變，而入口氣體的溫度Th、壓力Ph、空氣密度ρh等均會直接影響到空壓機出口氣壓Pcp的大小。因此，研究不同高度下的氣體特性對變高度下的空壓機建模十分重要。在海拔高度小于11 km時，Th及Ph與海拔高度h有如下關(guān)系[29]:

Th=T0+Lt(h-h0)

(1)

(2)

式中：T0=293 K;P0=101 325 Pa;h0=0 m;R為理想氣體常數(shù)；g0為重力加速度常數(shù)；Ma為空氣的摩爾質(zhì)量；Lt為溫度變化率。

空壓機入口空氣密度的大小直接影響到空壓機的能量轉(zhuǎn)換效率。水平面處的空氣密度約為1.29 kg/m3，海拔高度h處的空氣密度ρh可表示為

(3)

式中：Z為壓縮系數(shù)。

1.1.2 變高度下空壓機特性

與容積式空壓機不同，離心空壓機的輸出流量與壓力間存在強耦合關(guān)系。假設(shè)氣體壓縮為等熵壓縮，考慮能量傳遞對空壓機性能的影響，變高度下空壓機的工作特性可表示為

Δht(ωcp,Wcp)=Δhs(ωcp,Wcp)+Δhi(ωcp,Wcp)+

Δhf(ωcp,Wcp)+Δhother(ωcp,Wcp)

(4)

(5)

(6)

式中：Δht為焓的總增加量；Δhs為實際焓增量；Δhi、Δhf分別為氣體在葉片入口處和擴散器處因入射損失、摩擦損失引起的熵增；Δhother為在壓縮過程中由其他損失引起的熵增，這里假設(shè)為總熵增的8%;ωcp為轉(zhuǎn)軸的旋轉(zhuǎn)速度；Pcp、Tcp.out分別為空壓機出口的氣體壓力和溫度；ψ為空氣機壓比；cp為定壓比熱容；cv為定容比熱容；γ=cp/cv為比熱容比,約為1.4。

式(4)中空壓機出口流量Wcp可根據(jù)供應(yīng)歧管氣壓Psm與空壓機壓比ψ得到：

(7)

式中：A1為空壓機葉輪眼面積;Lc為空壓機氣體傳輸管道長度。

滯止焓的表達式為

(8)

(9)

壓縮過程中的摩擦損失與雷諾數(shù)Re相關(guān)，不同條件下雷諾數(shù)的變化對空壓機的輸出性能會產(chǎn)生較大影響。在變高度條件下，摩擦損失Δhf與Re的關(guān)系如式(10)～式(12)所示：

Δhf=Δhfi+Δhfd

(10)

(11)

(12)

式中：Δhfi和Δhfd分別為在葉輪和擴散器處的摩擦損失；li和ld分別為平均葉輪管道長度、平均擴散器管道長度；Di和Dd分別為葉輪處平均液壓管道直徑、擴散器處的平均液壓管道直徑；α2b為相對于切線方向的擴壓器進口角；A1為流通面積；f為平均雷諾數(shù)的摩擦系數(shù),

f=0.316 4(Re)-0.25

(13)

(14)

其中：D2為葉輪尖端直徑；u2為葉輪葉片出口的圓周速度；μ2為葉輪出口氣體溫度為T2時的動力黏度，假設(shè)為常數(shù)20.5×10-6kg·s/m2;ρ2為葉輪出口氣體密度。

基于式(1)～式(14)，給出了空壓機在0、2、3 km處的MAP圖，分別反映了在相應(yīng)高度下空壓機轉(zhuǎn)速、質(zhì)量流量以及出口氣壓間的動態(tài)關(guān)系，如圖2所示。其中，通過實驗測量獲得的空壓機地面數(shù)據(jù)在圖2(a)中給出。

由圖2可看出，所建立的空壓機特性模型與實驗測量數(shù)據(jù)有較好的吻合度，且隨著高度的增加，空壓機輸出氣體的質(zhì)量流量及壓力范圍逐漸變窄，空壓機正常工作范圍縮小。在海拔0、2、3 km 高度下，最大輸出質(zhì)量流量分別為0.024 2、0.019 5、0.017 5 kg/s，相應(yīng)的最大輸出壓力分別為1.732×105、1.391×105、1.251×105Pa。

圖2 不同高度下空壓機MAP圖

1.1.3 空壓機驅(qū)動電機建模

空壓機由每一時刻兩兩導(dǎo)通，三相六狀態(tài)運行的無刷直流電機(BLDCM)驅(qū)動，如圖3所示。

圖3 無刷直流電機原理圖

BLDCM三相電壓平衡方程為

(15)

式中：ua、ub、uc和ia、ib、ic以及ea、eb、ec分別是相電壓、相電流和各相的反電動勢；r=ra=rb=rc為定子電阻；L為定子電感；Vn為電機三相繞組中性點的電壓。采用星型連接，三相反電勢和三相電流間存在以下關(guān)系：

|ea|=|eb|=|ec|=Cφωcp

(16)

ia+ib+ic=0

(17)

式中：Cφ為磁鏈。

BLDCM三相理想反電動勢對應(yīng)的電流波形、霍爾信號如圖4所示。

圖4 理想反電動勢、電流及霍爾信號

BLDCM的動態(tài)特性方程為

(18)

(19)

(20)

式中：J為轉(zhuǎn)動慣量；τcm為無刷直流電機產(chǎn)生的電磁轉(zhuǎn)矩；τcp為空壓機的負(fù)載轉(zhuǎn)矩；ηcp為空壓機的效率。

1.2 歧管動態(tài)建模

1.2.1 供應(yīng)歧管建模

供應(yīng)歧管入口流量為空壓機的輸出流量Wcp，供應(yīng)歧管出口流量為Wsm.out，根據(jù)能量守恒和質(zhì)量連續(xù)方程，供應(yīng)歧管的壓力表達式為[30]

(21)

式中：Ra為空氣氣體常數(shù)；Tcp.out為空壓機出口氣體溫度；Vsm為供應(yīng)歧管體積；Tcp.out與空壓機壓比和效率有關(guān)，即

(22)

式中：ηcp為空壓機效率。

Wsm.out與Psm和陰極壓力Pca有關(guān)，即

Wsm.out=Ksm.out(Psm-Pca)

(23)

1.2.2 返回歧管建模

返回歧管的壓力與陰極出口流量Wca.out和返回歧管出口流量Wrm.out相關(guān)，表達式為

(24)

陰極出口流量的表達式為

Wca.out=Kca.out(Pca-Prm)

(25)

式中：Trm為返回歧管氣體的溫度，假設(shè)與燃料電池電堆溫度Tst相同；Vrm為返回歧管體積；Kca.out為返回歧管的流量系數(shù)。

返回歧管出口流量Wrm.out是關(guān)于返回歧管壓力Prm和末端環(huán)境氣壓Ph的函數(shù)：

(26)

式中：ADrm為返回歧管閥門的開口面積。

1.3 燃料電池電堆建模

1.3.1 陰極氣壓建模

陰極氣壓建模描述的是陰極內(nèi)部氣體的動態(tài)變化特性，其中，氧氣和氮氣是陰極流道中主要氣體成分，氧氣和氮氣分壓可根據(jù)氣體的質(zhì)量守恒和理想氣體定律得到[31]:

(27)

(28)

Pca=PO2+PN2+Psat

(29)

式中：Vca為陰極體積；MO2和MN2分別為氧氣和氮氣的摩爾質(zhì)量；WO2.in和WN2.in分別為陰極入口處氧氣和氮氣的質(zhì)量流量；WO2.out和WN2.out分別為陰極出口處氧氣和氮氣的質(zhì)量流量；WO2.react為陰極消耗氧氣的質(zhì)量流量；Psat為飽和水蒸汽壓力。

1.3.2 燃料電池電堆電壓建模

燃料電池電堆模型由N片燃料電池單體串聯(lián)組成，電堆輸出電壓為[32-34]

Vst=NVfc=N(Enernst-Vact-Vohm-Vcon)

(30)

式中：Vfc為燃料電池單體輸出電壓；Enernst為能斯特電動勢；Vact為活化損失過電勢；Vohm為歐姆損失過電勢；Vcon為濃差損失過電勢。能斯特電動勢以及3種損失過電勢的表達式為

Enernst=1.229-8.5×10-4×(Tfc-298.15)+

(31)

Vact=V0+Va(1-e-c1j)

(32)

(33)

Vohm=j·Rohm

(34)

式中：Tfc為燃料電池工作溫度，假設(shè)與燃料電池電堆溫度Tst相同；V0為電流密度為零時的壓降；Va與溫度和氧氣分壓有關(guān)；Rohm為歐姆內(nèi)阻；j和jmax分別為電流密度和最大電流密度；c1、c3均為常數(shù)；c2與飽和水蒸氣壓力和燃料電池工作溫度有關(guān)，電堆各參數(shù)的取值及表達式見附錄A。

2 控制器設(shè)計

PEMFC系統(tǒng)是一個具有非線性、強耦合性的多輸入多輸出系統(tǒng)，其系統(tǒng)效率和動態(tài)性能受到陰極供氣系統(tǒng)的直接影響。本文的控制目標(biāo)為基于所建PEMFC陰極供氣系統(tǒng)模型，避免氧饑餓，并保持PEMFC系統(tǒng)的高工作效率。為使控制效果最優(yōu)，本文提出了針對過氧比和陰極壓力的雙閉環(huán)控制，控制系統(tǒng)結(jié)構(gòu)如圖5所示。

如圖5所示，陰極壓力基于分?jǐn)?shù)階PIλDμ控制，通過調(diào)節(jié)返回歧管閥門開口面積來控制陰極壓力保持在給定值，過氧比控制由分?jǐn)?shù)階PIλDμ控制外環(huán)和MPC控制內(nèi)環(huán)組成，該控制策略簡寫為FFM。首字母F代表陰極壓力分?jǐn)?shù)階控制，后者F代表過氧比分?jǐn)?shù)階控制，M代表電機轉(zhuǎn)矩模型預(yù)測控制，并在圖5中分別給出了模型和控制模塊所對應(yīng)的公式。

圖5 FFM 控制系統(tǒng)框圖

2.1 陰極壓力與過氧比控制器設(shè)計

2.1.1 分?jǐn)?shù)階PIλDμ控制器

以陰極壓力分?jǐn)?shù)階PIλDμ控制器為例，結(jié)構(gòu)框圖如圖6所示。圖中，分?jǐn)?shù)階PIλDμ控制器的時域表達式為

圖6 陰極壓力控制結(jié)構(gòu)

u(t)=Kpe(t)+KiD-λe(t)+KdDμe(t)

(35)

式中：u(t)為返回歧管閥門開口面積；e(t)為陰極壓力參考值與系統(tǒng)輸出的差值ΔPca。

2.1.2 分?jǐn)?shù)階PIλDμ控制器設(shè)計

設(shè)計一個分?jǐn)?shù)階PIλDμ控制器需要確定5個參數(shù)：Kp、Ki、Kd、λ和μ。為使所設(shè)計的分?jǐn)?shù)階PIλDμ控制器控制效果最優(yōu)，在傳統(tǒng)參數(shù)整定方法——主導(dǎo)極值點法的基礎(chǔ)上，采用了一種新的參數(shù)整定方法——極點階數(shù)搜索法。并以陰極壓力控制為例，給出了基于極點階數(shù)搜索法的分?jǐn)?shù)階PIλDμ控制器設(shè)計的具體步驟，并以時間乘絕對誤差積分準(zhǔn)則(ITAE)作為控制參數(shù)整定的性能指標(biāo)。

步驟1估計比例增益Kp。

Kp的估計值可根據(jù)系統(tǒng)穩(wěn)態(tài)誤差ΔPca與理想給定值Pca.ref之比Et得到：

Kp=[100/(Et0Pca.ref)-1]a0

(36)

式中：a0為被控系統(tǒng)傳遞函數(shù)中的系數(shù)，防止Kp過大導(dǎo)致系統(tǒng)不穩(wěn)定。Et0和Et存在以下關(guān)系：

(37)

步驟2確定極點S1,2。

系統(tǒng)極點主要對系統(tǒng)的閉環(huán)響應(yīng)過程產(chǎn)生影響，當(dāng)閉環(huán)系統(tǒng)主要受到一對主導(dǎo)極點的影響時，該對主導(dǎo)極點用以下形式給出[35]

(38)

式中：σ為極點穩(wěn)定度;ξ為阻尼比。

當(dāng)步驟1和步驟2中的Kp和S1,2確定后，采用階數(shù)搜索方法來確定控制器中其余的4個參數(shù):Ki、Kd、λ和μ。

步驟3確定積分階次λ和微分階次μ的第i(1≤i≤3)次搜索范圍。

當(dāng)i=1時，μ∈[0,2.0],λ∈[0,1.5]，搜索步長Δ1=0.1。

當(dāng)i>1時，μ∈[max(Δi,μi-1-2Δi-1),μi-1+2Δi-1],λ∈[max(Δi,λi-1-2Δi-1),λi-1+2Δi-1]，搜索步長Δi=Δi-1/10。

步驟4確定i次搜索時的Kd、μ、Ki、λ。

將步驟1～步驟4中得到的Kp、S1,2和λ、μ的第i次搜索范圍中的值依次代入特征方程中，并令其等于零，即可解得每組的控制參數(shù)Kd、μ,Ki、λ。根據(jù)系統(tǒng)性能指標(biāo)(ITAE)來判定所得參數(shù)是否滿足要求，若有一組參數(shù)滿足指標(biāo)要求，則結(jié)束搜索；若均不滿足指標(biāo)要求，則記錄第i次搜索中最接近指標(biāo)要求的μi和λi，使i=i+1，若i≤3則返回步驟3，若i>3則進行步驟5。

步驟5重新搜索5個參數(shù)。

若系統(tǒng)階躍響應(yīng)時的ITAE較大或調(diào)節(jié)時間T較長，則適當(dāng)增大Kp，若系統(tǒng)階躍響應(yīng)的超調(diào)較大，則適當(dāng)減小Kp，并重復(fù)執(zhí)行步驟2～步驟4。具體的分?jǐn)?shù)階PIλDμ控制器設(shè)計流程如圖7所示。

圖7 最優(yōu)分?jǐn)?shù)階PIλDμ控制設(shè)計流程

文中選用ITAE優(yōu)化準(zhǔn)則的數(shù)學(xué)定義為

(39)

式中：r(t)為閉環(huán)系統(tǒng)給定參考值；y(t)為閉環(huán)系統(tǒng)實際輸出值。

2.2 電機驅(qū)動轉(zhuǎn)矩控制器設(shè)計

2.2.1 電流預(yù)測模型建立

當(dāng)無刷直流處于換相時刻時，a、b、c相電壓與直流母線電壓的關(guān)系為

(40)

δ值與非換相相、續(xù)流電流大小以及導(dǎo)通相開關(guān)管開通關(guān)斷有關(guān)，假設(shè)δx、δy、δz分別對應(yīng)非換相相、續(xù)流電流相和導(dǎo)通相。在電機運行6狀態(tài)下，非換相相δx、續(xù)流電流相δy和導(dǎo)通相δz因運行狀態(tài)的不同，均可表示電機a、b、c三相中的任意一相。

因此，狀態(tài)變量為相電流ia、ib、ic以及轉(zhuǎn)速ωcp，控制變量為開關(guān)狀態(tài)δa、δb、δc以及轉(zhuǎn)矩τcm。

上管調(diào)制：當(dāng)由c+b-轉(zhuǎn)換為a+b-拍時(即VT5、VT6導(dǎo)通→VT1、VT6導(dǎo)通)，

ea=-eb=E

(41)

1)ic≠0A(即c相電流續(xù)流)，若VT1導(dǎo)通，則δz=δa=1、δx=δb=0、δy=δc=0。若VT1關(guān)斷，則δz=δa=0、δx=δb=0、δy=δc=0。使用基爾霍夫定律得到當(dāng)VT1處于導(dǎo)通和關(guān)斷狀態(tài)時a、b、c相電壓方程為

VT1導(dǎo)通時:

(42)

VT1關(guān)斷時:

(43)

將式(15)、式(40)～式(43)離散化得到

(44)

(45)

式中：ib(k+1)、ib(k)分別為非換相相k+1時刻和k時刻電流值；Ts為采樣周期；on和off分別表示換相開關(guān)管開通及關(guān)斷狀態(tài)。

2)ic=0 A(即c相電流續(xù)流結(jié)束)，若VT1導(dǎo)通，則δz=δa=1、δx=δb=0、δy=δc=0。若VT1關(guān)斷，則δz=δa=0、δx=δb=0、δy=δc=0。使用基爾霍夫定律得到當(dāng)VT1處于導(dǎo)通和關(guān)斷狀態(tài)時a、b相電壓方程為

VT1導(dǎo)通時:

(46)

VT1關(guān)斷時:

(47)

將式(15)、式(41)、式(46)、式(47)離散化后得到

(48)

(49)

下管調(diào)制：當(dāng)由a+b-轉(zhuǎn)換為a+c-拍時(即VT1、VT6導(dǎo)通→VT1、VT2導(dǎo)通)，

ea=-ec=E

(50)

1)ib≠0A(即b相電流續(xù)流)，若VT2導(dǎo)通，則δx=δa=1、δy=δb=1、δz=δc=0。若VT2關(guān)斷，則δx=δa=1、δy=δb=1、δz=δc=1。使用基爾霍夫定律分別得到當(dāng)VT2處于導(dǎo)通和關(guān)斷狀態(tài)時a、b、c相電壓方程為

VT2導(dǎo)通時:

(51)

VT2關(guān)斷時:

(52)

將式(15)、式(40)、式(50)～式(52)離散化得到

(53)

(54)

2)ib=0A(即b相電流續(xù)流結(jié)束)，若VT2導(dǎo)通，則δx=δa=1、δy=δb=0、δz=δc=0。若VT2關(guān)斷，則δx=δa=1、δy=δb=0、δz=δc=1。使用基爾霍夫定律得到當(dāng)VT2處于導(dǎo)通和關(guān)斷狀態(tài)時a、c相電壓方程為

VT2導(dǎo)通時:

(55)

VT2關(guān)斷時:

(56)

將式(15)、式(50)、式(55)、式(56)離散化后得到：

(57)

(58)

基于包含逆變器數(shù)學(xué)模型的預(yù)測模型，有限集模型預(yù)測控制分別考慮了在導(dǎo)通相電流續(xù)流(iy≠0 A)和續(xù)流結(jié)束(iy=0 A)兩種情況下，導(dǎo)通相功率開關(guān)管分別開通或者關(guān)斷狀態(tài)時，電機非換相相未來時刻電流有4種可能預(yù)測值。同時，當(dāng)電機在ab、ac、bc、ba、ca、cb相導(dǎo)通的6種不同運行條件下，下一時刻預(yù)測值也是不同的，故每個周期共有24種有限情況。

2.2.2 滾動優(yōu)化

根據(jù)2.2.1節(jié)建立的電流預(yù)測模型，為減小非換相相電流波動，最小化非換相相電流以及對電流幅值進行約束，設(shè)計代價函數(shù)：

(59)

(60)

取min={J1,J2}，決定下一時刻導(dǎo)通相開關(guān)管的on或off狀態(tài)。

因無刷直流電機的極對數(shù)p=1，運行頻率為0～4 667 Hz。為保證模型與控制策略的可靠性，仿真速度應(yīng)大于50倍電機頻率，故文中設(shè)置模型預(yù)測的采樣頻率為250 kHz，步長Ts=4×10-6。文中采用的模型預(yù)測為一步預(yù)測，即預(yù)測時域為1。

圖8 模型預(yù)測控制系統(tǒng)框圖

3 仿真結(jié)果分析

為驗證建立的系統(tǒng)模型和提出的FFM控制方法，基于無人機飛行任務(wù)模擬過程，對控制的快速動態(tài)響應(yīng)特性和抑制轉(zhuǎn)矩脈動特性在MATLAB/Simulink上進行了仿真分析。同時，與PPM控制(首字母P代表陰極壓力PID控制、后者P代表過氧比PID控制、M代表電機轉(zhuǎn)矩MPC控制)和PPP控制(首字母P代表陰極壓力PID控制、中間P代表過氧比PID控制、后者P代表電機轉(zhuǎn)矩PID控制)進行了對比。燃料電池、空壓機及驅(qū)動電機和3種控制方法的詳細(xì)參數(shù)請見附錄A。

無人機飛行任務(wù)模擬過程如圖9所示。其中，陰極壓力與過氧比給定值的設(shè)置依據(jù)與不同高度下空壓機的工作特性有關(guān)，在20 s后，巡航負(fù)載電流小，空壓機輸入流量低，導(dǎo)致空壓機實際工作點在喘振線左側(cè)，為保證空壓機始終工作在正常工作區(qū)域，將過氧比參考值提高至5，同時將陰極壓力參考值降至地面氣壓1 atm(1 atm=101 325 Pa)。

圖9 無人機飛行任務(wù)模擬

3.1 陰極壓力閉環(huán)控制

3種控制下陰極壓力的控制效果如圖10所示，由圖可見，3種控制方法均可實現(xiàn)將陰極壓力控制在參考值附近，且穩(wěn)態(tài)誤差較小可忽略。在20 s參考值發(fā)生改變時，F(xiàn)FM調(diào)節(jié)時間相比PPM減小了0.5 s,僅為PPM調(diào)節(jié)時間的37%，PPM和PPP兩種控制下的陰極氣壓無明顯差異；在30 s和34 s時，因受空壓機流量控制的耦合影響，陰極氣壓分別有一個短暫減小和短暫增大的過程，但陰極氣壓減小和增大的幅度及響應(yīng)時間都較小，可忽略，因此圖中未給出這兩個時間的局部放大圖。

圖10 陰極氣壓控制效果

調(diào)節(jié)返回氣管閥門開口面積的大小是實現(xiàn)陰極氣壓控制的有效途徑之一，閥門開口面積越大，陰極氣壓越低。在20 s時，因陰極氣壓參考值減小，導(dǎo)致閥門開口面積突然增加，以期在短時間內(nèi)降低陰極氣壓跟蹤上給定值，閥門開口面積的變化情況如圖11所示。

3.2 過氧比閉環(huán)控制

3種控制方法下過氧比的控制效果如圖12所示。其中圖12(a)是過氧比的整體控制效果，可見在無人機飛行狀態(tài)發(fā)生改變時，3種控制方法均可將過氧比穩(wěn)定在給定值附近，F(xiàn)FM相比PPM和PPP，在調(diào)節(jié)時間上具有明顯的優(yōu)勢，而3種控制方法的超調(diào)量和穩(wěn)態(tài)誤差均差異不大。

圖12(b)為無人機飛行狀態(tài)由爬升到3 km 巡航時的局部放大圖，此時無人機需求功率急劇減小，且因3 km處空壓機特性的限制，需對過氧比做出相應(yīng)的改變。由圖可見，F(xiàn)FM控制的過氧比可以在較短的時間內(nèi)跟蹤上給定值，且相比PPM和PPP分別減少了75%和80%的調(diào)節(jié)時間，具體調(diào)節(jié)時間如表1所示。

圖12(c)為無人機飛行狀態(tài)由3 km巡航到下降時的局部放大圖，此時無人機需求功率繼續(xù)減小。由表1可知，F(xiàn)FM在調(diào)節(jié)時間上同樣具有最優(yōu)的控制效果，相比于PPM，其調(diào)節(jié)時間減小了60%，PPM和PPP的控制效果無明顯差異。

圖12(d)為無人機飛行狀態(tài)由下降到2 km 巡航時的局部放大圖，此時無人機需求功率增加，3種控制策略的調(diào)節(jié)時間如表1所示，F(xiàn)FM控制的調(diào)節(jié)時間僅為PPM控制的28%，為PPP控制的50%。

表1 3種控制方法過氧比控制的調(diào)節(jié)時間

過氧比定義為進入陰極的氧氣流量與反應(yīng)消耗氧氣流量的比值，因此，過氧比的控制與空壓機的流量緊密相關(guān)，圖13給出了空壓機流量的控制效果。

圖13 空壓機流量控制效果

電機驅(qū)動轉(zhuǎn)矩的脈動會導(dǎo)致其運行轉(zhuǎn)速的不穩(wěn)定，嚴(yán)重時更會造成電機的損壞。本文為減小轉(zhuǎn)矩脈動對電機的損害，提出了對轉(zhuǎn)矩的MPC控制?；贔FM、PPM和PPP這3種控制下的電機驅(qū)動轉(zhuǎn)矩如圖14所示，由圖可見，F(xiàn)FM和PPM控制下的電機轉(zhuǎn)矩特性無明顯區(qū)別，但相比PPP，轉(zhuǎn)矩脈動減小了接近50%，有效降低了轉(zhuǎn)矩脈動對驅(qū)動電機帶來的損害。

圖14 驅(qū)動轉(zhuǎn)矩控制

圖15為FFM、PPP和PPM控制下的a相電壓電流，對比可知，F(xiàn)FM和PPM控制下的a相電流脈動，僅為PPP控制的50%左右。因此，通過對圖14和圖15的對比分析可知，所提出的電機轉(zhuǎn)矩MPC控制算法可有效降低轉(zhuǎn)矩的脈動。

圖15 3種控制方法下的a相電壓、電流

文中無人機動力系統(tǒng)包含鋰電池和PEMFC兩種能源裝置，鋰電池僅用于啟動時空壓機的供電，PEMFC作為能源裝置為不同飛行狀態(tài)下的無人機提供需求功率，PEMFC系統(tǒng)輸出功率的大小與無人機飛行狀態(tài)緊密相關(guān)。圖16和圖17分別給出了3種控制策略下電堆電壓和功率的控制效果，由圖可見，無人機起飛及爬升階段、巡航階段和下降階段電堆輸出功率分別為9.4、4.4、3.5 kW，電堆電壓相應(yīng)為52、67、72 V，可滿足圖9中無人機在起飛和爬升階段的大功率需求。

圖16 電堆輸出電壓

圖17 電堆輸出功率控制效果

4 結(jié) 論

1)針對無人機用燃料電池在不同飛行高度下面臨的環(huán)境差異，考慮空壓機入口處氣體的溫度、壓力、密度以及雷諾數(shù)等影響，構(gòu)建了跨高度離心式空壓機模型、分析了不同高度下的空壓機特性MAP圖，在所建空壓機模型的基礎(chǔ)上完成了跨高度燃料電池陰極供氣系統(tǒng)建模及控制方法的研究。

2)基于無刷直流反電勢特征構(gòu)建了高速離心式空壓機驅(qū)動電機模型，同時針對電機驅(qū)動轉(zhuǎn)矩采用了有限集模型預(yù)測控制，并通過仿真驗證了MPC控制方法可有效降低驅(qū)動轉(zhuǎn)矩脈動。

3)所提出的FFM控制能夠根據(jù)陰極氣壓及過氧比的變化，及時調(diào)整陰極背壓閥開口面積和電機驅(qū)動轉(zhuǎn)矩，進而通過電機驅(qū)動轉(zhuǎn)矩調(diào)節(jié)空壓機轉(zhuǎn)速，以在短時間內(nèi)跟蹤上過氧比及陰極壓力參考量，具有可調(diào)參數(shù)范圍大、動態(tài)響應(yīng)速度快、電機轉(zhuǎn)矩脈動小等特點。所提出的FFM控制可實現(xiàn)燃料電池在0～3 km寬工況條件下運行，并滿足無人機不同飛行狀態(tài)下的功率需求。

附錄A

空壓機及驅(qū)動電機參數(shù)見表A1。

表A1 空壓機及驅(qū)動電機參數(shù)

Table A1 Parameters of air compressor and brushless motor

符號定義單位數(shù)值ψ地面最大壓力比1.7Wcp地面最大流量g/s24ηsmax地面最大等熵效率70m重量g600Pe額定功率W1 000p極對數(shù)1J轉(zhuǎn)動慣量kg·m25.5×10-7C?磁鏈Wb6.6×10-3Ie額定相電流A2.7r定子電阻Ω0.9L定子電感μH160σ滑移因子0.875A1流通面積m27.9×10-5β1b葉片入口角度(°)28β2b轉(zhuǎn)子葉片角度(°)80r1平均誘導(dǎo)半徑m0.010 3

燃料電池參數(shù)見表A2。

表A2 燃料電池參數(shù)

Table A2 Fuel cell parameters

符號定義單位數(shù)值N單電池個數(shù)100γ空氣比熱容比1.4Tst電池堆溫度K313.15Tatm環(huán)境溫度K298.15Vsm供應(yīng)岐管體積m30.01Vca陰極體積m30.01Vrm返回岐管體積m30.005Mv水蒸氣摩爾質(zhì)量kg/mol0.018MN2氮氣摩爾質(zhì)量kg/mol0.028MO2氧氣摩爾質(zhì)量kg/mol0.032Ma空氣摩爾質(zhì)量kg/mol0.029Φatm環(huán)境空氣的相對濕度0.5Ksm.out供應(yīng)管道孔口常數(shù)kg/(s·Pa)2×10-5Kca.ou陰極出口常數(shù)kg/(s·Pa)1×10-5mv.ca.max陰極側(cè)水蒸氣的最大質(zhì)量kg/mol0.0018yO2.atm氧氣摩爾分?jǐn)?shù)0.21yN2.atm氮氣摩爾分?jǐn)?shù)0.79

3種控制方法參數(shù)見表A3。

表A3 3種控制方法參數(shù)

Table A3 Parameters of three control methods

方法陰極壓力控制KpKiKdλμFFM0.0041×10-61×10-40.60.9PPM0.002 51×10-61×10-4PPP0.002 51×10-61×10-4方法過氧比控制KpKiKdλμFFM20.11×10-80.50.01PPM15200PPP25200方法電機轉(zhuǎn)矩控制KpKiKd預(yù)測時域仿真步長FFM14×10-6PPM14×10-6PPP0.90.851×10-4

燃料電池電堆參數(shù)為

V0=0.279-8.5×10-4(Tfc-298.15)+

Va=(-1.618×10-5Tfc+1.618×10-2)·

(-5.8×10-4Tfc+0.5736)

c1=10,c3=2,jmax=2.2 A/cm2

b11=0.051 39,b12=0.003 26,b2=350,

λm=14