基于偽衛星測試的并行計算災害救援中心選定模型構建

鞠森森，錢若霖

(陜西工業職業技術學院 土木工程學院，陜西 咸陽 712000)

0 引言

特殊的地質環境、不利的地形條件，導致我國各地區泥石流災害頻頻發生。當泥石流災害發生時，如何對應急資源進行科學、高效的調配，防止災害的進一步擴散，在理論和實踐上都有很大的作用。

泥石流災害救援中心的優化選址與救援時間滿意度是突發泥石流災害應急救援考慮的首要因素。其中，泥石流災害救援中心選定環節是最關鍵的一個環節，緊急物資能否準確、高效地運送到災區，直接影響到救災工作的成功和失敗。泥石流災害救援中心選定問題，即利用科學方法決定救援中心的地理方位，實現路徑最短的快速災害救援[1]。泥石流災害點面積廣、成災迅速，注重救援中心選址快速便利性，以便對其進行有效救援，成為災害修復工作的重點。

對此，有學者提出了相關研究。如，文獻[2]提出以覆蓋滿意度和經濟成本為目標的選址模型，構建自適應遺傳算法，并結合方向性提高搜索速度，通過自適應交叉與變異算子得到最優解。該方法不能根據道路被泥石流毀壞條件下進行再次尋優，適應度不足。文獻[3]提出以救援效率和成本為目標，構建多階段決策模型。分析每段路線的救援效率，構建多階段選址-路徑優化模型，利用改進快速非支配遺傳算法求解。該方法對混合交叉算子的限制不足，在求解時得到多個解，難以判定最優方案。

針對上述模型存在的問題，本文構建一種基于偽衛星測試的并行計算災害救援中心選定模型。使用基于并行計算的泥石流運動模擬方法，了解泥石流匯流運動過程，并將泥石流堆積作為泥石流全局運動，采用離散元算法建立泥石流堆積模型，將泥石流流團運動速率從三維變換為二維，預測泥石流災害圍線范圍；構建救援中心選定優化模型，在精度因子基礎上獲取救援選址的中心區域點，實現泥石流災害救援中心選定。

1 問題描述

在發生重大泥石流災害時，必須在最短的時間內，合理地選定災害救援中心，使其在有限的運輸條件下，盡快選擇合適的運輸方法，最大限度地滿足所需的物資運送時間和數量需求。泥石流災害救援中心選址問題是NP-Hard問題。可以用雙層模型來描述：上層模式是在最短的時間內，將救援物資從救援中心運送到對應的緊急轉運站，而不超出救援中心的距離和指定運輸站的的選址問題；下層模式是指在滿足一定的配送量要求條件下，以更短的時間、更低的運輸成本，將救援物資從緊急轉運站運送到不同的緊急情況點，而選擇的緊急轉運站則是兩種模式的連接點，由于上層模式的目標具有全局性，最優參數值依賴于自己的參數變量、約束條件，也要考慮下層的目標值，因此，必須及時向上層傳遞目標值；在此基礎上，上層決策對下層模式的路徑規劃產生了一定的影響，通過這種反饋機制，可以使上下兩層動態互動，從而達到兩層模式的綜合優化。

因此，在泥石流災害發生后，要根據救災物資的需求，既要考慮到救災物資的距離，又要綜合考慮各種物資的類型，以及所需物資的運輸費用，從而選擇一個合適的救援中心，并制定一個最佳的物資分配方案，以確保各個中心在最小的代價下，最大限度地滿足災民的需求，確保救災的效率，降低資源的消耗，從而得到系統的整體優化。問題描述如圖1所示。

圖1 災害救援中心選址-調度示意圖

由圖1可知，泥石流災害救援中心選址優化問題可以描述如下：

設有n個大概率發生災害的受災需求點，對于受災需求點k，其在預期內發生泥石流災害的概率為pk，在應急救援隊伍到達之前，災害擴散造成的損失函數為dk(t)，可容納的最大救援隊數量為qk，k=1，2，...，n，設有m個備用的救援中心選址，中心為O，運行成本為bo，救援階段反映函數為Tt，每次救援的成本為cl，l=1，2，...，m，從救援中心O到達受災需求點k所需的時間為tkl，k=1，2，...，n，l=1，2，...，m，為找出分布最好的救援中心的開設計劃和救援計劃，從而使預期損失和成本之和最小化，需要構建基于并行計算的泥石流運動學模型，確保泥石流救援中心的最優選址。

2 基于并行計算的泥石流運動學模型

本文在建立簡化泥石流堆積模型時，使用離散元+GPU算法實現泥石流堆積并行計算模型架構。運用離散元算法建立的泥石流運動公式，以顆粒之間的觸碰模型為主要研究重點。本文使用的觸碰模型是彈簧-阻尼器模型，將其法向接觸力描述為：

(1)

其中：fi，s是顆粒i與顆粒j撞擊的互斥力，顆粒在撞擊時，fi，d是顆粒阻尼力，k表示彈簧彈性指數，η是阻尼系數，d為顆粒直徑，vij是對應方位與兩個顆粒之間的對照速率，rij表示顆粒i與顆粒j之間的距離。

通過計算機模擬泥石流，利用地形特征對得到的運動數據資料，計算出泥石流的相關運動數據結果。結合流體力學模型法，并應用流團模式對泥石流的運動進行數值仿真。其中以堆積扇面為基礎，利用二維運動方程對其進行數值仿真。此模式具有簡單、高效的特點。在采用流團模式進行泥石流運動的計算機仿真時，必須對其進行時空離散，在此基礎上，將泥石流的運動過程劃分為若干個時間步長，將其空間劃分成相同的網格，利用統計學的方法，將各流團的運動參數表示為各個網格，從而得出泥石流的整體運動規律。將泥石流分為不同的多個粒子集合，這些粒子的體積和總的泥石流體積有很大的差異，因此可以認為是一大批顆粒在進行堆積扇的運動。在模型操作中，要將各個流團的運動速度從三維向二維轉換。

因此，在進行模型推算時，將流團空間速率V變換成二維速度方位內的速度[4-5]，使用迭代方法推算流團每個時段的速率，并將其進行簡化形變，得到流團迭代公式：

(2)

使用GPU迭代方法計算流團每個時段的速率，并將其進行簡化形變，得到如式(3)的流團迭代公式：

(3)

可以將式(3)直接應用在運算過程中，在時間域求解中，每一時間步長都需要進行一次迭代，以求出堆積過程中流團的泥石流速率[6-7]，摩阻坡降Sfx、Sfy與底面坡降Ssx、Ssy的運算過程為：

(4)

其中：τ、γm、μ、B、h、kcu表示泥石流參變量。在模擬泥石流運動形態的過程中，計算出流團在每個時間步長中的速率參變量[8-9]，同時也要獲取不同時段各個流團所處的方位參變量。已知第k個流團在n時段的方位坐標，那么迭代求解第k個顆粒在n及n+1時段的方位關系，具體表示為：

(5)

推算網格運動參量是一個相對動態過程，即獲取某個時段n的網格運動參數。在計算式，要將全部計算區域分割為正方形網格，并將其固定，然后按照運動流團某時段經過n的流團運動參數得到其網絡運動參量，運動參量分別是速率、泥深、動能及動量值[10-11]。

運用網絡內全部流團位于X、Y向速度的合速度值描述各個網格的泥石流速度值，將該過程記作：

(6)

通過式(6)計算流團位置，從而可以確定流團所在的柵格位置。在相同的流團網格中，利用流團體積，計算當前網格泥深，具體如式(7)所示：

hk=Δc*w/s

(7)

式中：hk表示在某一時刻第k個網格的泥深，w表示網格中流團數量，Δc表示單流團體積，s表示網格總面積。

在計算流圖動量和動能時，采用基于網格速度和網格體積的方法，由式(8)實現：

(8)

式中，Ek、Pk表示第k個網格在某一時刻的動能及動量，ρ表示泥石流密度，C表示網格內泥石流的體積，Vk表示第k個網格某一時刻的泥石流速度。

由前述求解過程可知，泥石流堆積模型中流團空間速率、泥石流速率被數值解耦，網格運動參量、流團位置求解、流圖動量和動能等計算密集型流程中沒有數據依賴關系的程序周期。因此，該模型的計算程序可以實現可并行性，實現多核心CPU和離散元的并行結構。當前，在 OpenMP基礎上實現多核 CPU的并行化處理已經比較成熟，所以在此不做詳細介紹。目前 GPU的并行處理技術有 CUDA、OpenCL、OpenACC等多種編程方式。在此基礎上，從軟件兼容和使用方便的角度出發，本文選擇了 OpenACC的并行編程方式，該模式是通過對編譯指令的識別和指令辨識，使其能夠產生并行的執行代碼，從而達到 CPU與離散元的異構并行運算。

并行計算與串行運算相關。該方法可以同時處理多條命令，以加快運算速度，并擴展問題的求解范圍，從而處理大量的、繁雜的運算。利用 GPU實現并行運算的方法有三大部分：①將所需要的存儲空間分布在主、設備兩個終端，并將其從主服務器拷貝到裝置；②通過設備端進行并行運算；③將計算的數據從設備上復制到主機，基于離散元+GPU算法是以 CPU為基礎的并行計算。

此外，如果將所有的泥石流都用離散元方法進行數值模擬，并通過CPU進行計算，該模型的規模將會非常有限。而在災害救援過程中，多個目標的集中救援和多個目標需求點，傳統方法很難得到精確的求解。因此，本研究在采用離散元法建立的泥石流運行累積模型中，所有的渲染運算都通過GPU完成，GPU則完成大量高密度的并行運算，CPU則承擔連續運算。通過 CPU/GPU的協作運算，可以對大規模泥石流匯流現象進行模擬，該研究結果可為面向流域泥石流匯流聚集風險區域的分析和救援中心選址提供更加精細的虛擬地理環境仿真。運用ThingJSglb三維模擬軟件模擬泥石流模型效果如圖2所示。

圖2 并行計算泥石流模型模擬結果

3 基于偽衛星測試的泥石流災害救援中心選定

在計算得到準確的泥石流模型后，將泥石流災害救援中心選定結合偽衛星測試技術，進行如下描述：假設偽衛星采集到具備m個可能的災害地點，關于災害地點j而言，其計劃期內產生災害的幾率是pj，在救援抵達之前，災害擴散擁有的損失函數是dj(t)，能夠容納的救援隊最多個數是Kj。若有n個備用的救援中心位置，中心i設立成本是Ci。從救援中心i抵達災害地點需要的時間是tij，決策人員希望尋找最佳分布的救援中心的設立方案和救援方法，使期望損失與成本為最低。允許多個救援中心聯合救援，一旦地點j發生災害，優先抵達該位置的救援中心次序是[1]，[2]，…，[Kj]，將聯合救援解析式定義為Uj(R[1](t)，R[2](t)，…，R[n](t))。如果災害地點j有泥石流現象發生，損失函數dj(t)會直線上升。在救援隊[1]在時段tj[1]抵達時，災害終止擴散，損失函數開始減小。救援隊[2]、[3]陸續抵達后，在多個救援小組的共同作用下，災害損失成快速下降趨勢[12-14]，因此可得到救援優化模型的模型參變量為：

i=1，2，…，n，j=1，2，…，m

(9)

(10)

將上述兩個公式使用以下數學模型進行描述：

(11)

(12)

(13)

xij，yi=0j=1，2，…，m

(14)

式中，在聯合救援函數U(°)內，如果x[k]j=0，那么證明R[k](t)不包括在U(°)中。

以上模型(9)～(12)是一個十分繁雜的0-1整數規劃模型。在目標函數式(9)內，第一項為救援中心的設立成本，針對以幾率pj產生的災害帶來的損失而言[15-17]，第一項表示救援成本，首個定積分為救援隊抵達前的損失，第二個積分為救援過程損失。約束式(10)代表參與救助的救援中心必須開設條件，式(11)為每個受災點所需的救援隊伍上限收斂，式(12)代表變量類別收斂。

假設某個時段用戶Ti探測到一個偽GPS衛星sj，同時構建以下測碼偽距觀察解析式：

(15)

(16)

不考慮大氣層的偏差影響，繼而獲取下列公式：

(17)

從上式可知，公式內含有四個未知參變量，因此最少構建四個相似的解析式[18-19]，用戶應同時觀察到四個偽衛星，滿足泥石流災害救援選址系統的精準定位。

矩陣Q內的矢量(li，mi，ni)為四顆定位偽衛星至測驗中心構成的向量和三個坐標軸夾角的方位余弦[20]。在系統內加入精度因子準則來評估選址定位結果，將平面精度因子與垂直精度因子描述為：

(18)

(19)

因為救援選址系統是以電子地圖作為解析目標，而電子地圖涵蓋的是地圖上各個點的二維位置，從上述的運算過程，方向余弦矩陣Q是推算精度因子的關鍵[21-22]，而方向余弦矩陣Q是通過參加定位的四顆偽衛星到用戶檢測點構成的四個向量和坐標系三個坐標軸夾角的余弦組成。

將檢測點O定位坐標系中心零點、四顆偽衛星及測量點組成四個向量OS1、OS2、OS3、OS4，將O當作球心組建一個單位球，獲得四個單位向量，且Z軸與地面保持垂直。

利用余弦定理求解OS2的水平方位角，OS2的方位角是其與OS1的水平夾角，記作：

(20)

通過上述計算，獲得原始測試點位置的HDOP值和VDOP值，繼而將其與評判標準進行對比，明確該位置是否是救援選址的核心區域點。

本文泥石流災害救援選址系統主要包含三個板塊，具體如圖3所示。

圖3 災害救援中心選定方法架構

衛星通信節點主要分布在災害地區附近，利用自組織模型構成通信網絡，以此進行數據收集與傳送。基站節點關鍵負責匯集節點傳送的有關數據，同時采用互聯網把數據輸送至服務器內。服務器處于救援中心，主要作用于信息儲存與加工，按照對應數據構建最佳的救援選址方案。

4 選址算法流程

由于災害擴散的損失遠高于救災的成本，因此距離災難現場最近的救援中心是參與聯合救援的救援中心的首選，其他參加救援的中心可以由離散元+GPU算法求出。采用離散元法與 GPU并行計算方法，建立了泥石流的流動堆積模式，以其堆積特征和災區為模擬參考，對泥石流運動進行模擬，根據模擬結果，預測泥石流覆蓋范圍，利用偽衛星技術計算災害擴散程度和易發點距離，將用戶選擇的位置與高度作為已知條件，以期望損失最低為約束，以精度因子為標準得到救援選址的核心區域點，可以避免陷入局部循環，最終可求得全局最優解。選址算法流程圖如圖4所示。

圖4 選址算法流程圖

5 仿真實驗

為了證明本文模型的可靠性，以某縣泥石流為例，對災害救援中心選定進行仿真分析。選取的泥石流流域面積5.66 km2，主溝長度為9 km，主溝比降142.5‰，最大高程為1 256 m，溝口高程為1 698 m，相對高度為1 242 m，流向由東向西。選址范圍在泥石流區域10公里半徑的圓內，設該地區有10個受災需求點。實驗平臺為Matlab 7.0仿真軟件，在處理器CPU為16 Core；Memory為128 GB；Disk為3TBBandwidth為1 000 Mb/s，內存為 DDR4 4 G的計算機上運行通過。

偽衛星測試得到的泥石流模擬效果如圖5所示。

圖5 泥石流仿真模擬效果圖

選址模擬圖結果如圖6所示。

圖6 選址模擬圖

通過圖4可以看出，黑色代表泥石流災害區域，淺黑色代表最優救援區域。本文給出不同的救援中心方案。通過不同方案均能快速到達最優救援區域。

為了驗證基于偽衛星測試的并行計算災害救援中心選定模型的定位準確性，將文獻[2]提出的基于自適應遺傳算法的應急物資儲備庫選址及物資調配優化模型、文獻[3]提出的帶模糊需求的多階段雙目標應急選址—路徑優化模型與本文模型進行對比，以救援中心選定定位精度為測試指標，三種模型救援中心選定定位精度對比示意圖如圖7所示。

圖7 救援中心選定定位精度對比

從圖7中可以看出，本文模型擁有較強的定位精度，隨著實驗次數的增加，其定位精度逐步穩定于90%以上，表明該模型具備較優的實用性，可以滿足災害后救援中心選定的及時選擇，為受災地區提供更快的救援效率。文獻[2]模型在實驗次數較少時，與本文模型的定位精度較為接近，但在實驗次數為5時，其選址精度僅有73%，而文獻[3]模型的救援選址定位精度最低，魯棒性較差。

為了驗證本文模型的有效性，為保證實驗數據的準確性，在同一試驗環境下進行試驗，確保實驗條件的一致性，以救援中心選址時間為測試指標，對本文模型、文獻[2]模型和文獻[3]模型的救援選址時間進行對比分析，對比結果如表1所示。

表1 救援中心選址時間對比

從表1中可知，文獻[2]模型的救援選址時長約為本文模型時長的2倍，而文獻[3]模型的救援選址時長約為本文模型時長的3倍，實用性均不高，無法提高災后救援的高效率執行。綜合判定，本文模型的救援選址效率優于傳統模型的救援選址效率，可廣泛應用于實際場景中。另外，根據上述結果可以看出，災害救援中心的選址在不同實驗次數下，由于救援地點的災情和需求的不同，呈現出動態變化的現象。在不同實驗次數下，根據災害救援中心的位置和救援地點的實際情況，合理選擇救援地點，最大限度地提高救援效率和救援地點的時間滿意度。

為了進一步驗證本文模型的可行性，采用本文提出的基于偽衛星測試的并行計算災害救援中心選定模型進行災害救援中心選址，并與文獻[2]模型和文獻[3]模型進行對比，3種模型獲得的災害救援中心位置坐標值如表2所示。

表2 災害救援中心選址結果

由表2中的災害救援中選址位置數值結果，繪制災害救援中心選定目標情況如圖8所示，其中圓形表示本文模型，三角形表示文獻[2]模型，實心方形表示文獻[3]模型，感知區域為核心區域，即災害救援中心的最優選址。

圖8 災害救援中心點最優解

對比三種模型的災害救援中心選址結果可知：感知區域的四個位置坐標分別為(3.53，3.11)、(5.76，3.11)、(3.53，8.09)、(5.76，8.09)，本文模型的5個災害救援中心選址點均在感知區域內部，且災害救援中心選址位置較近，方便受災人員轉移，本文模型的災害救援中心點最優解位置坐標為(4.23，5.12)；而文獻[2]模型和文獻[3]模型的災害救援中心選址點距離感知區域的最優位置坐標為(8.59，9.15)和(12.32，13.45)，且文獻[2]模型和文獻[3]模型的災害救援中心選址點均不在感知區域內部，其原因是未解決線性不收斂問題。經上述測試結果表明，本文模型的災害救援中心選址性能綜合最優。原因是利用偽衛星技術計算災害擴散程度和易發點距離。

為驗證本文構建一種基于偽衛星測試的并行計算災害救援中心選定模型中，離散元+GPU算法的收斂性，設置相同的最大迭代次數20，得到3種算法收斂曲線的對比結果如圖9所示。

由圖9可以看到本文提出的離散元+GPU算法，收斂速度優于其他兩種對比算法，上述實驗結果均表明本文算法的救援中心選址性能遠超于另外兩種對比算法。原因是該離散元+GPU算法運用離散元算法建立的泥石流運動公式，以顆粒之間的觸碰模型為核心要點，實現泥石流堆積并行計算模型架構。結合偽衛星測試技術，在系統內加入精度因子準則來評估選址定位結果，繼而將其與評判標準進行對比，明確該位置是否是救援選址的核心區域點，按照對應數據構建最佳的救援選址方案。同時，該算法采用增強了局部搜索能力，本文算法能迅速地將數據輸送至服務器內，進而表現出較好的優化性能。

結合上述內容，與其他常用于解決這兩類問題的兩種文獻模型做比較，得到的對比結果如表3所示。

表3 不同模型實驗性能的結果

從表3結果分析可以看出，文獻[2]模型和文獻[3]模型的迭代次數分別為20次和25次，均不能獲得全局最優解；而本文模型的迭代次數為15次，本文針對中心選址問題，該方法的并行運算采用了離散單元+ GPU的方法，能夠得到最優解，而且具有很高的效率。將兩種方法結合起來，可以得出一套科學的災害救援中心選址方案。

綜上所述，本文建立的模型具有極大的柔性，無論是救援地的時間滿意度，還是災害救援中心的選定，在實際泥石流災害救援中都是非常必要的，可以保證災害救援各個時期內最大限度地滿足救援地的物資需求，可以獲得全局最優解并且求解效率較高，同時使得救援地的時間滿意度和救援中心位置最優，為災害救援中心選定提供了科學依據。

6 結束語

為了增強災害救援的及時性，構建一種基于偽衛星測試的并行計算災害救援中心選定模型。通過分析固相泥石流顆粒之間的受力，采用離散元方法構建泥石流運動公式，將泥石流體分解為多流團顆粒集合，運用并行計算獲得泥石流運動模擬結果，準確評估泥石流災害具體情況；使用救援中心選定模型和偽衛星測試手段，完成并行計算災害救援中心選定模型構建，并增強救援選址模型整體性能，為真實場景下泥石流災害救援提供有效幫助。經實驗驗證可知，本文模型擁有較強的定位精度，隨著實驗次數的增加，其定位精度逐步穩定于90%以上，本文模型的5個災害救援中心選址點均在感知區域內部，且災害救援中心選址位置較近，方便受災人員轉移。因此，本論文的研究對緊急情況下的災害救援中心選定具有一定的現實意義。由于災害突發事件的影響因素較多，其優化模型的目標和約束條件相對簡單，沒有考慮到救援過程中的道路運輸條件、救援半徑、救援人員的應急反應速度、應急救援中心的周邊環境等，因此還有待于實踐的具體改進。