常 征，范瀚文，張聆曄

(大連海事大學(xué) 交通運(yùn)輸工程學(xué)院，遼寧 大連 116026)

0 引 言

隨著海上貿(mào)易的不斷發(fā)展，船舶交通流愈發(fā)復(fù)雜，船舶碰撞、擱淺、失火等突發(fā)事件屢見(jiàn)不鮮。以海上重大溢油事故為例，據(jù)ITOPF(international tanker owners pollution federation)2019年油輪溢油統(tǒng)計(jì)分析報(bào)告指出，雖然全球海上溢油事故數(shù)呈現(xiàn)逐年下降的趨勢(shì)，但2010年—2018年總計(jì)發(fā)生溢油量超過(guò)7 t事故數(shù)超60余起。因此，在事故發(fā)生后，如何科學(xué)的響應(yīng)應(yīng)急需求，快速的送抵應(yīng)急物資是應(yīng)急救援的重中之重。為此，亟需解決在一定資源的約束下，實(shí)現(xiàn)對(duì)應(yīng)急物資儲(chǔ)備庫(kù)的選址優(yōu)化。

應(yīng)急物資儲(chǔ)備庫(kù)選址問(wèn)題屬于設(shè)施選址問(wèn)題的研究范疇，自1911年Weber提出單一倉(cāng)庫(kù)設(shè)施選址問(wèn)題以來(lái)，國(guó)內(nèi)外學(xué)者對(duì)設(shè)施選址問(wèn)題展開(kāi)了大量的研究，取得了豐碩的成果。叢雯婧等[1]綜合考慮了選址決策中的風(fēng)險(xiǎn)、效率以及經(jīng)濟(jì)性的目標(biāo)，求解了臺(tái)風(fēng)情境下的應(yīng)急選址問(wèn)題；張玲等[2]考慮了應(yīng)急救援在不同階段的需求變動(dòng)并建立了選址-分配模型，利用情景松弛的迭代算法進(jìn)行求解，證明了算法和模型的有效性；以上研究大多將選址模型應(yīng)用在陸上選址中，也有相關(guān)學(xué)者針對(duì)海上應(yīng)急物資儲(chǔ)備庫(kù)的建設(shè)問(wèn)題進(jìn)行研究。Y.F.AI等[3]考慮了政企聯(lián)合存儲(chǔ)的海上物資應(yīng)急儲(chǔ)備庫(kù)選址-分配模型，并設(shè)計(jì)了兩階段啟發(fā)式算法進(jìn)行求解，通過(guò)算例證明了政企聯(lián)合存儲(chǔ)是穩(wěn)定的雙贏策略；張聆曄等[4]考慮了應(yīng)急資源調(diào)度中及時(shí)性的要求，建立了風(fēng)險(xiǎn)不確定下海上應(yīng)急物資儲(chǔ)備庫(kù)選址魯棒優(yōu)化模型。

多重能力(多等級(jí))設(shè)施選址問(wèn)題是設(shè)施選址問(wèn)題的重要分支。國(guó)外相關(guān)學(xué)者已將多重能力等級(jí)選址模型應(yīng)用在醫(yī)療救助、通訊網(wǎng)絡(luò)、火災(zāi)防治、制造生產(chǎn)選址等方面[5-6]。K.HONORA等[7]將多重能力等級(jí)選址模型應(yīng)用在發(fā)展中國(guó)家的貧困農(nóng)村醫(yī)療設(shè)施選址問(wèn)題中，在模型中考慮了時(shí)效性和公平性的要求，通過(guò)實(shí)例驗(yàn)證了模型和算法的有效性；J.MICHAEL等[8]從災(zāi)難發(fā)生后群眾需求的差異性角度出發(fā)，設(shè)計(jì)了多重能力等級(jí)P-median模型，結(jié)合GIS系統(tǒng)以及空間優(yōu)化布局理論，采用CPLEX方法進(jìn)行模型求解；常征等[9]構(gòu)造了帶有多重能力等級(jí)的內(nèi)陸港選址模型，并設(shè)計(jì)遺傳算法進(jìn)行求解，結(jié)果表明在考慮多重建設(shè)能力下的選址方案有效節(jié)約了投資成本；陸相林等[10]運(yùn)用改進(jìn)的多重能力選址模型優(yōu)化了消防設(shè)施配置，對(duì)北京市房山區(qū)展開(kāi)實(shí)證研究驗(yàn)證了模型和算法的有效性。

綜上所述，以往針對(duì)于應(yīng)急物資儲(chǔ)備庫(kù)選址問(wèn)題的研究大多在固定備選點(diǎn)建設(shè)等級(jí)下展開(kāi)，多應(yīng)用于陸上設(shè)施選址問(wèn)題。設(shè)施建設(shè)是一項(xiàng)系統(tǒng)工程，具有投資回收期長(zhǎng)、投資金額巨大的特點(diǎn)。同時(shí)，建設(shè)海上應(yīng)急物資儲(chǔ)備庫(kù)的選址設(shè)施模型中需要考慮高風(fēng)險(xiǎn)的水域的分布情況[11]。當(dāng)應(yīng)急需求產(chǎn)生時(shí)，不僅需要考慮應(yīng)急物資儲(chǔ)備庫(kù)與應(yīng)急需求點(diǎn)的空間距離，同時(shí)也要考慮應(yīng)急物資儲(chǔ)備庫(kù)的服務(wù)能力。為了更準(zhǔn)確地對(duì)應(yīng)急需求點(diǎn)與應(yīng)急物資儲(chǔ)備庫(kù)的互動(dòng)關(guān)系、服務(wù)能力、空間距離等因素進(jìn)行綜合考量，筆者擬引入引力模型對(duì)傳統(tǒng)選址模型進(jìn)行改進(jìn)。構(gòu)建多等級(jí)應(yīng)急物資儲(chǔ)備庫(kù)選址-分配模型，并設(shè)計(jì)改進(jìn)蝙蝠算法進(jìn)行求解。

1 問(wèn)題描述

假設(shè)在某一海上事故高發(fā)的海域附近建立若干應(yīng)急物資儲(chǔ)備庫(kù)，以滿足應(yīng)急救援需要。

備選應(yīng)急物資儲(chǔ)備庫(kù)地理位置已知，具備不同的可選建設(shè)等級(jí)。建設(shè)級(jí)別越高，應(yīng)急救援服務(wù)質(zhì)量越好。若應(yīng)急需求點(diǎn)與應(yīng)急物資儲(chǔ)備庫(kù)的距離超過(guò)該建設(shè)等級(jí)下的最大服務(wù)半徑，則該應(yīng)急物資儲(chǔ)備庫(kù)無(wú)法針對(duì)該需求點(diǎn)提供有效的應(yīng)急救援服務(wù)。

假定應(yīng)急需求點(diǎn)的地理位置、應(yīng)急物資需求、風(fēng)險(xiǎn)等級(jí)通過(guò)歷史統(tǒng)計(jì)推理可獲得。問(wèn)題的目標(biāo)是如何在有限的資金以及備選點(diǎn)的限制下確定建設(shè)應(yīng)急物資儲(chǔ)備庫(kù)的位置以及等級(jí)，滿足應(yīng)急需求點(diǎn)的應(yīng)急物資儲(chǔ)備需求。

2 選址模型建立

2.1 引力模型

引力模型起源于牛頓萬(wàn)有引力定律，即任意兩個(gè)質(zhì)點(diǎn)之間的引力與空間距離的大小成反比，與質(zhì)量大小成正比。引力模型作為經(jīng)典模型的通用規(guī)律，已經(jīng)在空間轉(zhuǎn)移、貿(mào)易流通、人文地理、場(chǎng)站選址等方面取得了較為成熟的應(yīng)用[12]。

考慮到海上事故的不可預(yù)知性、成因復(fù)雜性，海上應(yīng)急救援往往要綜合考慮應(yīng)急救援力量與需求點(diǎn)之間的相互“吸引關(guān)系”。即應(yīng)急物資儲(chǔ)備庫(kù)的建設(shè)等級(jí)越高(容量大小、應(yīng)急救援船只數(shù)目、應(yīng)急救援人員業(yè)務(wù)水平等)越高，則其對(duì)應(yīng)急需求的服務(wù)水平越好。除此之外，海上應(yīng)急物資儲(chǔ)備庫(kù)對(duì)應(yīng)急需求點(diǎn)的吸引程度也存在“空間距離遞遠(yuǎn)，吸引能力遞減”的現(xiàn)象。

故筆者引入牛頓式基本引力模型，刻畫(huà)應(yīng)急物資儲(chǔ)備庫(kù)與應(yīng)急需求點(diǎn)的“吸引力”為：

(1)

式中：Mi為應(yīng)急物資儲(chǔ)備庫(kù)的“質(zhì)量”，用于表征應(yīng)急物資儲(chǔ)備庫(kù)的建設(shè)水平，即應(yīng)急物資儲(chǔ)備庫(kù)的服務(wù)能力；Nj為需求點(diǎn)的“質(zhì)量”，用于表征需求點(diǎn)的風(fēng)險(xiǎn)等級(jí)，即應(yīng)急需求點(diǎn)的風(fēng)險(xiǎn)程度；dij為應(yīng)急物資儲(chǔ)備庫(kù)與應(yīng)急需求點(diǎn)的歐式距離；λ為萬(wàn)有引力常數(shù)，為了便于衡量應(yīng)急物資儲(chǔ)備庫(kù)與應(yīng)急需求點(diǎn)的“吸引力”，筆者取值為1；β為距離衰減系數(shù)，因筆者采用歐氏距離刻畫(huà)兩點(diǎn)之間的空間距離，故取值為1。

2.2 模型建立

建立海上多等級(jí)應(yīng)急物資儲(chǔ)備庫(kù)選址-分配模型為：

(2)

(3)

(4)

(5)

(6)

(7)

(8)

(9)

(10)

(11)

模型中式(2)和式(3)為目標(biāo)函數(shù)，式(2)為最大化系統(tǒng)引力；式(3)為最小化系統(tǒng)費(fèi)用，包括應(yīng)急物資儲(chǔ)備庫(kù)的建設(shè)費(fèi)用與應(yīng)急物資存儲(chǔ)過(guò)程中產(chǎn)生的存儲(chǔ)費(fèi)用；式(4)為所有應(yīng)急物資儲(chǔ)備庫(kù)的建設(shè)成本約束；式(5)為所有應(yīng)急需求點(diǎn)至少被覆蓋一次，即應(yīng)急需求點(diǎn)的全覆蓋；式(6)為建設(shè)應(yīng)急物資儲(chǔ)備庫(kù)不能超過(guò)規(guī)定的個(gè)數(shù)；式(7)為每個(gè)備選應(yīng)急物資儲(chǔ)備庫(kù)有多個(gè)建設(shè)等級(jí)，但最多只能選擇建設(shè)一種；式(8)為應(yīng)急物資儲(chǔ)備量不能超過(guò)k建設(shè)等級(jí)下的最高存儲(chǔ)量限制；式(9)為滿足應(yīng)急需求點(diǎn)的應(yīng)急物資需求；式(10)～式(11)為決策變量。

模型所涉及的參數(shù)以及決策變量如下所述：

1)參數(shù)集合：I為應(yīng)急物資儲(chǔ)備庫(kù)集合；J為應(yīng)急需求點(diǎn)集合；Dj為應(yīng)急需求點(diǎn)的應(yīng)急物資需求量；Tk為建設(shè)應(yīng)急物資儲(chǔ)備庫(kù)建設(shè)規(guī)模等級(jí)集合k=1,2,3；Cik為應(yīng)急物資儲(chǔ)備庫(kù)i在k建設(shè)等級(jí)下的建設(shè)成本；fik為應(yīng)急物資儲(chǔ)備庫(kù)i在k建設(shè)等級(jí)下的單位存儲(chǔ)成本；B為建設(shè)應(yīng)急物資儲(chǔ)備庫(kù)的最大數(shù)量；Hj={I|dij≤rik}為應(yīng)急物資儲(chǔ)備庫(kù)i在k建設(shè)等級(jí)下能夠覆蓋的應(yīng)急需求點(diǎn)J集合；P為投資額限制；Rik為應(yīng)急物資儲(chǔ)備庫(kù)i在k建設(shè)等級(jí)下的最大服務(wù)半徑。

3 算法設(shè)計(jì)

蝙蝠算法是一種利用模擬自然界蝙蝠捕食過(guò)程的新型仿生智能算法。蝙蝠捕食的原理為：向環(huán)境發(fā)射超聲波，當(dāng)遇到獵物時(shí)，接受反射回來(lái)的回聲，利用雙耳效應(yīng)以及回聲的強(qiáng)度變化分析獵物的具體方位、大小、移動(dòng)方向等，進(jìn)而迅速的捕捉獵物。蝙蝠發(fā)出的恒定頻率根據(jù)搜尋的目標(biāo)大小而改變，一般持續(xù)時(shí)間為8～10 ms，頻率為25～150 kHz，響度為110 dB。

蝙蝠算法具有強(qiáng)收斂性、搜索范圍大、具有潛在分布式等優(yōu)點(diǎn)，因此被廣泛應(yīng)用在多目標(biāo)優(yōu)化的問(wèn)題中。求解的模型是NP-Hard問(wèn)題。為了彌補(bǔ)蝙蝠算法初始種群隨機(jī)生成，不具備覆蓋整個(gè)解空間的缺陷，選用全局混沌蝙蝠優(yōu)化算法對(duì)問(wèn)題進(jìn)行求解。

3.1 基于混沌映射的種群初始化

蝙蝠的個(gè)體以迭代的方式進(jìn)行更新，種群內(nèi)部的全局最優(yōu)解是一個(gè)未知量，初始化種群的規(guī)模、范圍、影響對(duì)后續(xù)求解的收斂速度、最優(yōu)解的質(zhì)量有巨大影響，因此初始化種群的覆蓋范圍很重要。有效的改進(jìn)方法就是全面覆蓋整個(gè)初始空間，對(duì)求解最優(yōu)解的質(zhì)量有推動(dòng)意義。筆者采用混沌映射的方法對(duì)種群進(jìn)行初始化操作，提升初始解空間的利用率。目前應(yīng)用較為廣泛的為logistic方法[13]，計(jì)算方法為:

(12)

(13)

式中：lw和uw為變量取值范圍的最大值和最小值。

3.2 自適應(yīng)網(wǎng)格密度估計(jì)

首先建立q維目標(biāo)空間的自適應(yīng)網(wǎng)絡(luò)，將目標(biāo)空間劃分為D1×D2×…×Dq的網(wǎng)絡(luò)，其中第w個(gè)目標(biāo)的網(wǎng)絡(luò)模長(zhǎng)計(jì)算方式為：

(14)

外部檔案成員的數(shù)量會(huì)隨著外部存檔中非占有解的個(gè)數(shù)變化而變化，因此目標(biāo)空間的網(wǎng)格數(shù)量將會(huì)隨之自適應(yīng)變化。非占優(yōu)解在自適應(yīng)的空間位置為：

(15)

式中：mod(·)函數(shù)為比值后的整數(shù)部分，根據(jù)式(15)可以計(jì)算出每個(gè)網(wǎng)格內(nèi)部包含的非占優(yōu)解的數(shù)量，將一個(gè)網(wǎng)格內(nèi)的非占優(yōu)解的數(shù)量作為粒子的分布密度，粒子分布以及對(duì)應(yīng)的分布密度如圖1。

圖1 種群密度估計(jì)Fig. 1 Population density estimation

3.3 全局最優(yōu)解的選擇策略

全局最優(yōu)解是一組Pareto最優(yōu)的解集，如果沒(méi)有采取有效的選擇策略，筆者采取自適應(yīng)的網(wǎng)絡(luò)密度和輪盤(pán)賭法進(jìn)行全局最優(yōu)解的選擇，可以使每一代的全局最優(yōu)解不完全相同，具體的算法為：

Step 1利用自適應(yīng)網(wǎng)格法評(píng)估粒子的個(gè)體密度。

Step 2將個(gè)體的密度等級(jí)為一致的粒子統(tǒng)一放到集合中。

Step 3將個(gè)體的密度等級(jí)除以所有的個(gè)體密度等級(jí)之和，采用輪盤(pán)賭法來(lái)選擇個(gè)體的密度等級(jí)。

Step 4在相應(yīng)的個(gè)體密度等級(jí)集合中隨機(jī)選擇粒子作為全局的最優(yōu)解。

4 算例仿真

假設(shè)在某事故高發(fā)海域附近建立應(yīng)急物資儲(chǔ)備庫(kù)以滿足該海域的應(yīng)急救援需求。為有效實(shí)施選址優(yōu)化，將算例中連續(xù)的水域單元離散化處理，抽象為互不重合的點(diǎn)。同時(shí)，采取一組正交的網(wǎng)格線將算例的求解區(qū)域進(jìn)行分割。為驗(yàn)證算法的有效性，筆者在50×50平面上，隨機(jī)生成8個(gè)擬定建立應(yīng)急物資儲(chǔ)備庫(kù)的備選點(diǎn)、25個(gè)需求點(diǎn)以及對(duì)應(yīng)的應(yīng)急物資需求，見(jiàn)表1和表2。假定備選應(yīng)急物資儲(chǔ)備庫(kù)在不同等級(jí)下的建設(shè)成本、單位存儲(chǔ)成本、服務(wù)能力見(jiàn)表3。

表1 應(yīng)急物資儲(chǔ)備庫(kù)備選點(diǎn)數(shù)據(jù)Table 1 Optional point data of emergency material reserve

表2 應(yīng)急需求點(diǎn)數(shù)據(jù)Table 2 Emergency demand point data

表3 多重等級(jí)應(yīng)急物資儲(chǔ)備庫(kù)建設(shè)參數(shù)Table 3 Multi-level emergency material reserve constructionparameters

采用MATLAB軟件進(jìn)行程序設(shè)計(jì)，初始化參數(shù)設(shè)置：蝙蝠數(shù)量為40，響度為0.1，脈沖發(fā)射率為0.85，最低頻率為0.1，最高頻率為0.85，最大迭代次數(shù)為200。計(jì)算得出118個(gè)解，考慮解的分布特點(diǎn)，筆者選出22個(gè)非劣解(表4)組成的帕累托前沿見(jiàn)圖2。

圖2 帕累托前沿Fig. 2 Pareto frontier

表4 改進(jìn)蝙蝠算法求解帕累托解集Table 4 Pareto solution set achieved by improved bat algorithm

由圖2可見(jiàn)，求解得到的帕累托前沿分布較為均勻，可以為決策者提供有效的決策參考。同時(shí)，系統(tǒng)費(fèi)用與引力值存在正相關(guān)關(guān)系，即隨著等級(jí)建設(shè)能力、應(yīng)急物資儲(chǔ)備量的增加，引力值隨之增長(zhǎng)。決策者可以根據(jù)自身偏好或者現(xiàn)實(shí)需求，選擇滿意解作為應(yīng)急物資儲(chǔ)備庫(kù)的選址方案。為了保證選址結(jié)果的相對(duì)客觀，筆者提出引力損失率和系統(tǒng)費(fèi)用增長(zhǎng)率指標(biāo)，以輔助決策者針對(duì)帕累托前沿進(jìn)行選擇。分別用LF和LC來(lái)表示方案的引力損失率和系統(tǒng)費(fèi)用增長(zhǎng)率，計(jì)算公式為：

(16)

(17)

式中：max(LF)和min(LF)分別為系統(tǒng)引力值的最大值和最小值；max(LC)和min(LC)分別為系統(tǒng)費(fèi)用的最大值和最小值；決策者可以根據(jù)對(duì)系統(tǒng)引力損失率以及對(duì)系統(tǒng)費(fèi)用增長(zhǎng)率的權(quán)衡來(lái)選擇最終的選址方案。

決策者可以根據(jù)引力損失率以及系統(tǒng)費(fèi)用增長(zhǎng)率的偏好程度選擇1個(gè)滿意的求解方案。若決策者可以接受20%的引力損失，那么決策者需要在引力損失小于等于20%的方案中選擇系統(tǒng)費(fèi)用增長(zhǎng)率最高的方案作為本身的滿意解。其他分析方案同理。筆者僅給出了當(dāng)引力損失率和系統(tǒng)費(fèi)用增長(zhǎng)率為20%、40%、60%的方案選擇結(jié)果作為決策參考，結(jié)果如圖2。

為了客觀分析模型的實(shí)用性，筆者求解了固定設(shè)施等級(jí)下的選址結(jié)果。決策者在考慮最優(yōu)解時(shí)，其實(shí)是一個(gè)多屬性決策問(wèn)題。為了便于比較，筆者采用逼近理想解法[14](technique for order preference by similarity to ideal solution，TOPSIS)計(jì)算備選方案與正理想方案、負(fù)理想方案的距離來(lái)判斷最優(yōu)方案，進(jìn)而選出問(wèn)題的最優(yōu)解。比較結(jié)果見(jiàn)表5。

表5 固定等級(jí)下選址結(jié)果Table 5 Site selection results under a fixed level

結(jié)果分析：在固定最小建設(shè)等級(jí)下，備選應(yīng)急物資儲(chǔ)備庫(kù)無(wú)法實(shí)現(xiàn)對(duì)應(yīng)急需求的全面覆蓋，不滿足應(yīng)急救援的實(shí)際需求，故引力值較??；在建設(shè)中等等級(jí)下，備選應(yīng)急物資儲(chǔ)備庫(kù)無(wú)法實(shí)現(xiàn)對(duì)應(yīng)急需求的有效覆蓋，即部分應(yīng)急需求點(diǎn)沒(méi)有被全部滿足應(yīng)急物資的儲(chǔ)備需求；在固定最高建設(shè)等級(jí)下，備選應(yīng)急物資儲(chǔ)備庫(kù)實(shí)現(xiàn)了對(duì)應(yīng)急需求的完全滿足，但是耗費(fèi)了大量的建設(shè)資金。相較于多重能力選址模型而言，引力值增加了8.49%，系統(tǒng)成本降低了19.94%。

綜上所述，考慮了備選點(diǎn)多重能力等級(jí)的選址模型可以提升應(yīng)急救援的效率，節(jié)省經(jīng)濟(jì)成本。最終確認(rèn)多重等級(jí)設(shè)施選址求解結(jié)果見(jiàn)表6。

表6 應(yīng)急物資儲(chǔ)備庫(kù)備選點(diǎn)與應(yīng)急需求點(diǎn)分配Table 6 Location selection of emergency material reserve andallocation of emergency demand points

5 結(jié) 語(yǔ)

筆者研究了多重能力等級(jí)海上應(yīng)急物資儲(chǔ)備庫(kù)選址問(wèn)題。以最大化系統(tǒng)引力值和最小化系統(tǒng)成本作為目標(biāo)，將同時(shí)反映應(yīng)急救援能力和空間距離特征的引力模型結(jié)合全覆蓋選址思想，建立了多目標(biāo)混合整數(shù)規(guī)劃模型，并設(shè)計(jì)改進(jìn)蝙蝠算法取得帕累托解集。結(jié)果表明：系統(tǒng)引力值與系統(tǒng)費(fèi)用之間呈現(xiàn)正相關(guān)的關(guān)系，通過(guò)引入引力損失率和系統(tǒng)費(fèi)用增長(zhǎng)率兩個(gè)參數(shù)，可以為決策者提供不同偏好下的決策。

為了驗(yàn)證多重能力等級(jí)選址模型的優(yōu)越性，筆者對(duì)比了固定能力等級(jí)下的選址結(jié)果。通過(guò)算例仿真驗(yàn)證在滿足應(yīng)急物資需求的前提下，多重能力等級(jí)選址模型可以有效的降低經(jīng)濟(jì)成本，提高應(yīng)急物資儲(chǔ)備庫(kù)的運(yùn)行效率。該方法可以為海上應(yīng)急物資儲(chǔ)備庫(kù)選址提供一定的參考。