楊新湦, 屈琮博, 王梓旭
(1.中國民航大學中歐航空工程師學院, 天津 300300; 2.中國民航大學空中交通管理學院, 天津 300300)
航空運輸系統是一個復雜的巨系統,航班不正常是每個航空公司不可避免會面臨的問題。航空公司大面積不正常航班恢復是困難問題,一般而言,航空公司不正常航班恢復可分為分階段恢復與一體化恢復。分階段恢是將整個恢復問題分為飛機恢復、機組恢復、旅客恢復三個階段逐步求解;一體化恢復則是綜合考慮飛機、機組、旅客因素。
對不正常航班恢復的研究多數集中于分階段恢復研究。Teodorovic等[1-4]在1984~1995年先后對飛機停運、機場關閉航班恢復問題、旅客恢復和機組恢復問題進行了研究。分別使用分支定界、動態規劃、啟發式算法進行了求解,奠定了不正常航班分階段恢復的基礎。吳剛等[5]提出一種改進列生成算法,每次迭代過程中加入多個列,并對加入的多個列需要滿足的條件進行了分析,最后給出的算例驗證了該方法的正確性和有效性。田倩南等[6]改進時空網絡算法,給出占優準則,有效減少了可恢復航線的組合數量,提升計算效率。

已有文獻中多采用時空網絡進行建模,時空網絡中不僅有描述航班銜接的航班邊,還有機場節點和時間節點。當進行一體化建模時導致模型過于復雜。Sherali等[12]系統介紹了航班網絡的使用,相比時空網絡模型復雜度低,但是不能將時間、地點等屬性表現在圖中,諸如飛機維修約束、機組執勤時間約束、客票取消操作不易表達。對于巡航速度控制的策略,Aktürk等[13]在飛機恢復問題中首次使用,并分析了巡航速度控制的優勢,Arkan等[14]基于此研究進一步將旅客行程加入模型,但是都沒有達到真正的一體化恢復。
基于以上文獻中的不足,現做出如下工作:①改進航班網絡的弊端,并設計一種基于廣度優先搜索的網絡生成算法;②建立巡航速度控制下的一體化恢復模型;③使用雙曲不等式修正模型,符合二階錐規劃條件后進行求解并通過算例說明巡航速度控制在一體化恢復中的作用。
針對航班網絡相較于時空網絡不能表達出飛機維修、調機、客票取消、加機組等操作的劣勢進行改進。傳統航班網絡如圖1所示,其中包括三類節點和三類弧,改進航班網絡中有四類節點和五類弧,新加入必經節點和必經弧、添加弧,具體如下。

sr為初始節點;tr為沉沒節點;Mr為必經節點集合;f、g為航班
源節點:表示飛機、機組或旅客初始狀態的點,有時間、地點屬性。地點屬性是飛機、機組、旅客行程的首發機場,時間屬性是首發時間。
末節點:表示飛機、機組或旅客的終止狀態的點,有時間、地點屬性。地點屬性是對飛機停場過夜機場的限制,時間屬性是過夜機場的宵禁時間。
航班節點:表示飛機、機組執行的航班或旅客在其行程中要乘坐的航班。
必經節點:表示飛機在特定時段在特定機場維修的節點,或機組規定在特定機場休息。
出發弧:連接源節點與航班節點的弧,表示飛機、機組、旅客執行首個航班。
沉沒?。哼B接航班節點與末節點的弧,表示飛機、機組、旅客執行完最后一個航班。
航班?。汉桨喙濣c之間的弧,表示飛機、機組、旅客執行完上游航班后執行下游航班。
必經?。哼B接航班節點與必經節點之間的弧,表示飛機維修、機組休息等。
添加?。罕硎竞桨嗷謴偷奶厥獠僮?,如加機組、調機操作(航班節點之間、源節點與末節點之間、源節點與航班節點、航班節點與末節點之間);乘客客票取消操作(源節點與末節點之間、航班節點與末節點之間)。
一體化恢復包含飛機、機組、旅客三種子網絡,用r表示。算法設計中,F為航班集合;closetimedown為機場開始關閉時間;closetimeup為結束關閉時間;curfew為機場宵禁時間;Ntemp為航班暫時存儲器中符合條件的航班f集合;Nnext為航班暫時存儲器中符合條件的航班g集合;Orig為航班g的始發地,Desf為航班f的目的地;ctfg為航班銜接時間;Nr為節點集合;Vr為弧集合;Er為添加弧集合。航班網絡生成算法由Generatesubnetwork(r),Generatepath(Ntemp)和Insert(Ntemp) 三部分組成。具體步驟如下:
Step1 初始化節點集合Nr,弧集合Vr,Ntemp=sr。
Step2 Generatesubnetwork(r),航班網絡生成主程序如下。
Procedure Generatesubnetwork(r)
Gr=(Nr,Vr)←Generatepath(Ntemp) #為第二部分程序
Vr=Vr∪Er
Return Gr=(Nr,Vr)
end Procedure
Step3 Generatepath(Ntemp),基于廣度優先搜索可銜接航班,即航班網絡路徑中的節點。
Procedure Generatepath(Ntemp)
WhileNtemp≠?
f←first element ofNtemp
ifclosetimedown≤[dtf,atf]≤closetimeuporatf≤curfew
Updatedtfandatf#修改航班的OD時間
ifdtf-sdtf> 最長延誤時間
continue
else
Nnext←(g∈FOrig=Desfandsdtg≥atf+ctfg)
#當前航班的可銜接航班
for eachginNnext
Ntemp=Ntemp∪{g}
ifg∈Nr
Insert(Ntemp) #第三部分節點/弧生成程序
else
if Desf=Desr
end if
end if
end for
end if
end if
end procedure
Step4 Insert(Ntemp),得到航班網絡中的所有弧與節點。
Procedure Insert(Ntemp)
以實施“美麗鄉村小康水”行動計劃為主抓手,貴州省全面推進民生水利建設,讓更多群眾在水利發展中得到實惠。2013年全省預計解決300萬農村人口及學校師生的飲水安全問題,實施中小河流項目129個,治理水土流失面積2 200km2,實施病險水庫治理271座,新增小水電裝機20萬kW,新增小型農田水利重點縣22個,新增有效灌溉面積73萬畝,全面完成年初預定的各項目標。
Nr=Nr∪Ntemp∪Mr
fori=1:len(Ntemp)
Vr=Vr∪{fi,fi+1}
end for
end procedure
Step5 對Vr進行分類得到出發弧、沉沒弧、航班弧、必經弧、添加?。粚r進行分類得到初始節點、航班節點、沉沒節點。
Step6 如果是機組子網絡,額外刪除不滿足機組執勤時間和飛行時間的弧與節點。
飛機在最大航程速度(maximum ragne cruise,MRC)下最省油,在巡航速度大于MRC速度時燃油消耗隨速度增加而增加。速度控制綜合考慮加速與減速研究速度控制對航班恢復的影響,減速考慮飛機的反常操作,按照文獻[15-16]的建議巡航速度調整范圍在MRC速度的[-10%,+10%]。將巡航速度納入模型需要速度與燃油消耗之間的關系以計算恢復成本,歐洲交通管理組織的飛機數據基礎(base of aircraft data, BADA)項目中開發的巡航階段燃料流模型,提出給定質量和高度,可以計算出飛機巡航階段的燃油消耗率fcr(v)(kg / min)與速度v(km / min)的關系[17]為

(1)
式(1)中:c1、c2、c3、c4分別為與飛機阻力、油耗系數、給定高度的空氣密度和重力加速度相關的系數,這些系數可以從399種飛機的BADA用戶手冊中獲得。給定燃油消耗率可以計算出巡航階段總燃油消耗,假設飛行距離固定為d,那么可以得到巡航總耗油為

(2)
飛機以MRC巡航最省油,但是實際運行中MRC速度接近反常操縱區,常采用經濟巡航或LRC巡航,假設航班的計劃飛行速度都為經濟巡航速度vECON(按照MRC巡航速度的1.02倍計算),那么巡航速度調整后,航班燃油成本ΔFuelcost變化為
ΔFuelcost=pfuel[F(v)-F(vECON)]
(3)
式(3)中:pfuel為單位燃油價格;F(vECON)為經濟巡航速度下總耗油。
飛機巡航階段不僅是飛行時間、燃油消耗最多的階段,也是污染物排放最多的階段。國際民航組織(International Civil Aviation Organization,ICAO)根據發動機廠家審定的數據,定義了標準起飛著陸循環(land take-off,LTO),包括起飛、滑行、爬升、著陸4個階段[18],在基礎排放模型中給出4個階段的基礎排放數據見表1。Tamara使用ICAO數據結合波音公司的BM2方法對英國一天的流量樣本的燃油消耗與二氧化碳排放進行了分析,并與英國最廣泛使用的NETCEN估算進行了對比,確定了其實用性[19]。從表1可以看出,CO2的排放指數與飛行階段無關,魏志強等[20]研究指出CO2的排放量與燃油消耗成正比,占飛機排放物總量的99%以上,同時也是節能減排關注的重點。巡航速度的改變使燃油消耗變化導致二氧化碳排放的因素在模型中也將予以考慮,二氧化碳排放量Δcarboncost可以按照式(4)計算。

表1 CFM56-7B26 發動機基準排放數據[18]
Δcarboncost=pCO2k[F(v)-F(vECON)]
(4)
式(4)中:pCO2為單位燃油二氧化碳排放量;k為二氧化碳指數。
在模型提出之前做如下假設。
(1)飛機除座位數之外的性能差異不大,即不同飛機對巡航速度調整的限制相同,航班被不同飛機執行若沒有速度調整其計劃巡航時間保持不變。
(2)飛機的計劃巡航速度設為MRC的1.02倍,且不考慮航線上風和溫度變化的影響,速度調整范圍在[-10%,+10%]。
(3)對于多段行程的旅客不能使其到達中間地點后無后續航班,在恢復期內對旅客行程受到擾動的旅客進行轉到其他航班調整,不受擾動的航班上的旅客不可調整后受擾。
(4)飛行中的非巡航部分如滑行、起飛、爬升、下降的總時間為40 min。
模型算法中的集合參數變量定義如下:
集合
ac飛機集合MidN中間節點集合
cr機組集合MN必經節點集合
ps旅客行程集合arcs弧集合
Rr∈ac,cr,pssrarcs出發弧集合
N節點集合trarcs沉沒弧集合
SN首節點sr集合midarcs中間弧集合
TN末節點tr集合cntarcs聯程航班弧集合
F航班集合psnumber安排到航班的人數之和
參數

atf航班f計劃到達時間tf航班f飛行時間
DTf航班f實際起飛時間Np行程p的旅客人數
ATf航班f實際到達時間Nf原航班上旅客人數
df航班f的巡航飛行距離S飛機可用座位數
ctfg中轉銜接時間etr飛機或機組的最早可用時間

客客票取消成本
Cfd航班f取消成本CCO2碳排放成本
δtf航班f巡航縮短時間Cfuel燃油成本
Cps旅客單位延誤成本
變量



對于每種網絡要滿足網絡的流平衡約束和節點閉合約束。流平衡約束指每個網絡的源節點到中間節點的弧有且只有一個,每個網絡的中間節點到沉默節點的弧有且只有一個,中間節點的進入弧與出發弧數量相等。節點閉合約束指若弧fg被執行則ac與cr中的弧fg必有且只有一個被執行。如約束(5)~約束(10)。約束(5)~約束(7)是經典圖論模型的基本約束,約束(8)表示若弧fg被經過則飛機、機組網絡中必同時有弧fg被執行,約束(9)與約束(10)表示航班f只能被一架飛機和一個機組執行,由于旅客行程可以將旅客分配到任一或多個OD相同的航班因此無此約束。

(5)

(6)

(7)

(8)

(9)

(10)
對于不同的飛機有飛機特性約束,如飛機座位數約束,旅客行程被安排到航班f時的人數之和不能大于執行這個航班飛機的座位數如約束(11)。被分配旅客行程的人數與旅客行程取消人數之和等于原旅客行程人數如約束(12)。不同的飛機對于速度的調整限度相同,但是由于各飛機的MRC速度不一樣因此巡航時間調整并不相同,不同航線上的速度調整限制不同,如約束(13)。約束(14)表示機組僅可以執行執照允許的飛機任務,ac(cr)表示可以執行飛機ac的機組集合,約束(15)表示飛機不可以執行規定之外的航線(如高原航線由特定飛機執飛),ac(f)表示飛機ac可以執行的航線集合。

(11)

(12)

(13)

(14)

(15)


(16)

(17)

(18)
0≤DTf-dtf≤480
tfd=ATf-atf,f∈F
(19)
0≤atf≤curfewf,f∈F
(20)
根據第2.1節的闡述,將速度變化帶來的燃油消耗變化納入約束中,為了減少模型中的變量將速度用距離除以巡航時間表示,系數c1、c2、c3、c4也相應改變,如約束(21)~約束(24),其中Tfc是計劃巡航時間,當不考慮巡航速度時,則沒有約束(21)~約束(24)。

(21)

(22)

(23)

(24)
模型的目標函數以最小化成本表示,成本分為航空公司延誤運營成本(Cop,包括航班運營延誤成本與機組空閑成本)、航班取消成本(Cfd)、旅客延誤成本(Cpd)、客票取消成本(Cpc)、速度變化帶來的額外燃油(Cfuel)與碳排放(Cco2)成本。民航局在《航班延誤經濟補償指導意見》中的賠償建議,延誤8 h以上給出的賠償標準已經接近于取消航班的賠償,所以我們認為航班取消成本是航班延誤8 h的運營成本;旅客客票取消成本是旅客延誤8 h經濟損失加上客票價格。

(25)
二階錐規劃(second-order cone programming,SOCP)可以視作線性規劃的推廣, 本質上是一種凸規劃, 解的最優性和計算高效性都有優良特性。利用現有的CPLEX等算法包可以獲得較好的求解結果, 諸多基于SOCP的優化問題可以在多項式時間內完成。本節用到相關定義如下。
(2)上鏡圖:函數f:Rn→R的圖像定義為{[x,f(x)]|x∈domf},它是Rn+1空間上的一個子集,則epif={(x,t)|x∈domf,f(x)≤t}是該函數上鏡圖[21]。還有定理,函數f是凸函數當且僅當epif是凸集。
(3)雙曲不等式:如u2≤v1v2,易證明雙曲不等式可寫成二階錐不等式‖(2u,v1-v2)‖≤v1+v2。
一體化恢復模型中約束(21)是非線性約束,是非線性不連續的,因此其上鏡圖是非凸的,約束非凸,模型中的目標函數與其他約束條件也受約束(21)的影響。為了簡化說明刪除了變量索引,約束(21)的上鏡圖可表示為
t≥df(c1λ1+c2λ2+c3λ3+c4λ4)
(26)
y4≤λ1t2y
(27)
y2≤λ2t
(28)
t2≤λ3y
(29)
t4≤λ4y2t
(30)


(31)
式(28)、式(29)是雙曲不等式,式(27)與式(30)都可以寫成如下雙曲不等式:
y2≤st,s2≤λ1y
(32)
t2≤sy,s2≤λ4t
(33)
由雙曲不等式可以寫為二階錐不等式的性質可知,模型可以重構為一個混合整數二階錐規劃問題,模型有線性目標函數,二階錐約束,重構后的模型可以使用CPLEX求解得到結果。
采用某航空公司一天的運營數據,數據包含120 個航班,26 架飛機,6 種機型,32 個機場,16 473名旅客。機型參數如表2所示。針對既有長時間機場關閉又有飛機故障設置兩個場景。設場景一是浦東機場從該日應開放時刻關閉至12:00,機尾號為

表2 機型的參數
5113飛機故障4 h,分別使用本文模型與不考慮巡航速度控制的模型(即去除速度相關的約束與目標成本)進行求解構成算例1、算例2。場景二是浦東機場從該日應開放時刻關閉至下午16:00,機尾號為5113飛機故障8 h,分別使用本文模型與不考慮巡航速度控制的模型(即去除速度相關的約束與目標成本)進行求解構成算例3、算例4。規定航班最長延誤時間為8 h;機組最長執勤時間為14 h,最長飛行時間為8 h;航班之間最短銜接時間為40 min;航班單位運營延誤按機型尾流強度分為,重型機4 167元/h,中型機2 916元/h,輕型機208元/h,本文算例中包含中型機 24架,重型機2 架;旅客單位延誤經濟損失,國內旅客50 元/h,國際旅客100元 /h;機組空閑成本50 元/min。算例信息如表3所示,求解結果如表4所示。

表4 求解結果

表3 場景算例信息
算例1、算例3與算例2、算例4對比可以看到巡航速度控制的最優成本小于不考慮速度控制。但當恢復場景復雜,受到擾動航班較多,考慮巡航速度控制,恢復成本減少更加明顯。圖2與圖3表示4個算例的旅客延誤分布與恢復成本分布,考慮巡航速度控制后算例2較算例1在0~1 h的延誤人數減少448人,旅客延誤成本減少20 690元,額外燃油成本6 224元,碳排放成本915元。算例4較算例3少取消1班航班,旅客行程取消人數減少52人,客票取消成本46 956元,旅客延誤成本減少8 643元。

圖2 旅客延誤分布Fig.2 Passenger delay distribution

圖3 成本分布Fig.3 Cost distribution
各機場的延誤情況如圖4與圖5所示。算例2較算例1將湛江與成都的延誤人數清零,在武漢機場的延誤人數從788 人降到637 人,武漢機場累計延誤時間下降95 min,浦東機場累計延誤時間下降22 min。算例4較算例3將成都、曼谷的延誤人數清零,武漢機場累計延誤時間下降57 min。從成本分布來看,巡航速度控制將等時間的旅客延誤成本轉化為燃油消耗與碳排放成本,兩者之間的差值使總成本減小,而這個差值的大小當受擾航班較多時會比較明顯。

圖4 各機場延誤人數Fig.4 Delays by airport

圖5 各機場累計延誤時間Fig.5 Cumulative delays at airports
以航班擾動恢復為目的,考慮飛機、機組、旅客一體化恢復。根據提出的航班網絡生成算法,應用改進后的航班網絡,將巡航速度控制策略應用于航班恢復建立一體化恢復建立模型。模型易于理解算法容易計算機編程實現,通過擾動場景模擬,可在短時間內獲得調整方案。結果表明巡航速度控制策略通過將等時間的旅客延誤成本轉化為燃油消耗與碳排放成本,在恢復場景復雜,受到擾動航班較多時可明顯減少恢復成本,對航班延誤恢復問題由一定實際意義。