王文濤,甘旭升*,吳亞榮,任謹慎,王翠香
(空軍工程大學,a. 空管領航學院;b. 基礎部,陜西 西安 710051)
隨著航空技術與信息技術日趨成熟,無人機系統技術發展迅猛,同時,憑借其低成本、高機動、操作便捷等特點,應用日益多樣化[1-2]。而低空空域(特別是400~500 m 高度以下)作為“低、慢、小”無人機運行場所,能夠滿足無人機多樣化應用(偵察監視、物流快遞、救援救災、農業耕種等)需求[3-6]。但低空無人機數目多、體積小、執行任務種類多樣,現有的空中交通管理技術與手段難以適用于低空目標。因此,如何保證在低空復雜態勢下,實現無人機安全和高效的運行,已成為全球工業界和學術界的共同挑戰[7]。
目前,為解決低空無人機安全高效運行問題,各國都普遍在無人機運行服務系統(UAV service system,USS)系統和飛行情報管理系統(flight information management system,FIMS)支撐下,構建無人機空中交通管理(unmanned aircraft system traffic management,UTM)體系,并將 UTM 作為管控無人機的主要手段。美國聯邦航空管理局(federal aviation administration,FAA)通過 FIMS,在 USS、ATM 和國家空域系統(national airspace system,NAS)之間的協調機制,并通過UTM 維持空中態勢感知,維持空域和交通運行的監管與運營,同時給空域用戶提供使用權限[8]。
風險評估主要是指利用科學的計算方法和模型,對潛在的或實際事故造成的影響后果進行定量估計,這是決定無人飛行系統能否安全融入國家空域系統的關鍵環節。為提升UTM 概念中的實時安全評估能力,推進無人機系統逐步融入NAS,近年來已經有一些學者對無人機低空運行風險進行評估。Giulio Avanzini 等[9]提出了一種適用于在人口稠密地區的無人航空系統風險評估程序。該方法根據無人機系統可靠性、人口密度以及撞擊點的位置來評估發生故障后對地撞擊風險。Stefano Primatesta 等[10]利用風險圖對無人機運行風險進行描述,該方法可以對無人機運行環境的風險程度進行量化。其中風險圖通過概率方法生成,包括人口密度、遮蔽因子、禁飛區和障礙物等層級,無人機每個運行位置都有其對應的量化風險值。P. Gon?alves等[11]提出一種基于Petri 網無人機安全評估模型,通過Petri 網顯示已識別的無人機頻率,以及無人機在出現故障后的反應能力,對無人機運行中操作授權過程進行優化。Achim Washington 等[12]基于貝葉斯方法,充分考慮系統安全評估過程中固有的不確定性,對風險評估與決策框架進行定義,并以無人機安全評估系統為研究案例,驗證了所提方法的合理性。Xuejun Zhang 等[13]考慮并分析了無人機地面撞擊操作和空中碰撞的2 個關鍵危害,通過計算滿足不同無人機類別的安全目標水平所需的系統可靠性,來構建多因素影響下的風險模型,從而實現地面安全風險評估。REN X H 等[14]為評估低空運行的無人機的空域風險,定義了無人機在城市運行中的第三方風險——無人系統對地面無關人員造成的風險。并且明確了風險的來源和對象。考慮到無人機墜毀、噪聲、機載攝像頭和地面環境因素,構建了無人機城市物流風險指數評價模型。
從當前研究現狀來看,低空無人機運行風險評估研究已有一些進展,但往往忽略復雜低空環境參數具有不確定性特征,以確定的數值進行仿真計算,往往與實際運行情況相差過大。基于以上考慮,本文提出一種考慮無人機墜落不確定性的風險評估方法。根據無人機下降過程中的受力模型,推導出無人機的墜落區域,并考慮無人機運行過程中的不確定性,計算出無人機墜落時的事故危害程度。
與有人機相比,無人機在性能參數、任務場景等方面都有顯著差異,因此在評估無人機墜落地面撞擊事件造成的危害程度時,需要結合UAV 自身特點,對事故發生危害性進行建模評估。如圖1 所示,從UAV 系統來看,與地面撞擊事件危害程度相關的因素包括無人機自身屬性、墜落影響面積等。從任務場景角度看,無人機運行區域內地面人口分布、屬性與遮蔽參數是地面撞擊事故需要重點考慮的因素。

圖1 UAV 地面撞擊事故安全風險要素Fig.1 Safety risk elements of UAV ground impact accident
當前,主要用UAV 地面撞擊事故死亡人數P 作為衡量事故等級的標準[15],表示為

式中:P( f ) 為被撞擊人員死亡率;Nexp為暴露在UAV 地面撞擊事故中的地面人數,可由UAV 任務場景下的人口分布密度ρ 與撞擊事故影響區域面積Aexp進行計算:

無人機系統由于運行故障,墜落形式通常為拋體運動,在下落過程中受多因素影響,其運動方程[16]為

式中:m,G 分別為無人機質量與重力;x,z 分別為下落過程中的水平與垂直方向的距離;t 為下落時間;Dx,Dz分別為無人機下落過程中所受的空氣阻力。

式中:CD為空氣阻力系數;ρa為空氣密度;g 為重力加速度;Ax,Az分別為無人機側向與垂向迎風面積。

式中:vx0,vz0分別為無人機水平和垂直方向上的初始速度;x0,z0分別為無人機水平和垂直方向上的初始位置。
無人機墜落時間為

式中:h0為無人機墜落時的運行高度。
將式(6)代入式(5)中,可得無人機在水平初速度v0和運行高度h0墜落條件下的墜落撞擊點位置。實際無人機運行過程中的高度精度和水平精度分別記為 δh0,δxz。
無人機撞擊致死率表示為

式中:ps為遮蔽參數,表示無人機運行環境下人群的暴露程度,取值范圍為(0,1],參數取值越大,同等撞擊條件下致死率越低,各類型區域遮蔽參數設置如表1 所示。α 定義為ps= 0.5 時,致死率達到50%需要的撞擊能量;β 定義為ps= 0 時,發生死亡所需的撞擊能量臨界值。Eimp指撞擊事故中無人機動中無人機水平方向與垂直方向的碰撞速度vx,vy可由式(5),(6)確定,表示為


表1 遮蔽參數Table 1 Masking parameters
由于無人機低空運行區域內人口分布、地面遮蔽效果、飛行高度等參數均不相同,每個運行位置點存在的風險值也不一致,根據式(1)~(5),可對三維空間內運行風險進行量化,每一個地理坐標上的位置都對應一個風險值[17]。構建 n × m × l 的三維風險矩陣R,其中,第k(k = 1,2,…,l)層高度所代表的二維風險矩陣為

式 中 :r(pi,j,k)(i = 1,2,…,n;j = 1,2,…,m;k =1,2,…,l)表示三維低空空域內位置處于 pi,j,k時的風險值。根據上文分析,可用UAV 地面撞擊事故死亡人數P 來代替。
無人機在低空運行空域所受多類型不確定因素影響,例如,各類導航、制導與控制(guidance,navigation,and control,GNS)誤差,地面撞擊事故死亡人數P 同樣也難以用確定的數值進行描述,特別是無人機墜落實際撞擊點將隨機分布在標準撞擊點周圍,其概率密度函數(probability density function,PDF)依賴于撞擊誤差的統計分布。因此地面撞擊區域的統計影響區域(statistical impact footprint,SIF)可通過 Monte-Carlo 方法來確定,圖 2 展示了無人機以水平初速度v0和高度h0飛行時,拋體墜落運動地面撞擊點二維示意圖。

圖2 無人機拋體墜落運動地面撞擊點二維示意圖Fig.2 Two-dimensional schematic diagram of the ground impact point of a UAV projectile falling motion
三維空域內無人機墜落撞擊點分布也有類似的規律,通過Monte-Carlo 模擬可以確定雙變量概率密度函數:

式中:Δx = x - μx,Δy = y - μy,其中 μx,μy分別為軌跡方向與垂直軌跡方向的位置平均值;ρxy為x 與y之間的相關系數;σx> 0,σy> 0 分別為軌跡方向與垂直軌跡方向的位置方差。
無人機墜落撞擊點Monte-Carlo 模擬分布如圖3所示,墜地撞擊點三維分布圖如圖4 所示。其中紅色曲線表示墜落撞擊點分別在x 軸和y 軸的分布,是立體三維分布在二維面上的投影。

圖3 無人機墜落撞擊點Monte-Carlo 模擬分布Fig.3 Monte-Carlo simulation distribution of the impact point of a UAV fall

圖4 墜地撞擊點三維分布圖Fig.4 Three-dimensional distribution map of the impact points of a fall
不同概率密度函數下的撞擊影響點分布與影響區域面積如圖5 所示,按照撞擊區域3σ 準則,為充分考慮無人機對地面不同區域的撞擊概率不同,準確衡量無人機在空域內某個運行位置的風險值,將撞擊點分布圖區分為 Sσ、S2σ與 S3σ3 種類別,若無人機墜落影響區域經過柵格化后,分為n 個基本計算單位,則風險值表示為

圖5 不同撞擊概率下影響面積示意圖Fig.5 Schematic diagram of affected area under different impact probabilities

式中:a,b,c 分別為 Sσ、S2σ- Sσ與 S3σ- S2σ3 種區域面積內覆蓋的柵格個數,滿足a + b + c = n。
低空空域環境內地面障礙物類型復雜多樣、人口密度高,且處于動態變化之中。選取108°54′44″N~108°55'13″N,34°14'28″E~34°14'41″E 作為無人機運行區域,如圖6 所示。該地區包含建筑、空地、綠化等多類型地形,同時人口密集,是典型的人口聚集區,適合用來衡量無人機運行風險。考慮到實際運行評估過程中的可操作性,確保無人機在高于地面障礙物的高度上運行,其不干擾民航飛機的運行,選擇該區域上空30 m 空域作為無人機運行高度進行風險評估。

圖6 無人機運行區域Fig.6 UAV operating area
選取大疆經緯M210-RTK 型無人機作為風險評估對象,該型號無人機的參數如表2 所示。

表2 大疆經緯M210-RTK 型無人機參數Table 2 Parameters of DJI Jingwei M210-RTK UAV
人口密度是衡量地面人員分布的重要指標,是無人機運行風險評估中的一個關鍵變量。由于區域內每日人口空間分布具有較強的規律性與動態性,因此本研究選取 2020-04-20T12:00:00 固定時刻的城市熱力圖,作為衡量人口分布的依據,結果如圖7 所示,將運行區域柵格化,形成80 × 41 個網格,網格中不同顏色代表不同的人口分布密度。

圖7 人口密度分布圖Fig.7 Population density distribution map
遮蔽參數是衡量地面不同設施對人員的保護效果的重要指標,如圖8 所示,運行區域柵格化,形成80 × 41 個網格,根據表1 中定義,網格中不同顏色代表不同的遮蔽參數數值。

圖8 遮蔽參數分布圖Fig.8 Distribution of shadowing parameters
如圖9 所示,A 點為無人機任務起始點,經過A-B-C-D 路線最終返回A 點,運行高度設置在所有地面設施高度之上,為簡化評估流程,假設無人機在任務全程保持運行速度為v0= 20 m/s,運行高度h0= 30 m 不變。風險評估所需的參數如表3 所示。

圖9 無人機任務場景Fig.9 UAV mission scenario

表3 所需參數取值Table 3 Values of required parameters
為衡量無人機運行過程中的危害性,需要對離散化的運行區域進行逐一分析,對于所選取的區域,規模為:740 m×401 m,則每個柵格網絡的規模大約為:10 m×10 m。現在對初始點狀態的無人機進行風險評估,其余路徑點類似。
初始點無人機速度為v0= 20 m/s,高度h0= 30 m,以墜地撞擊點縱向分布為例,根據統計信息,其均值服從均值為42.51 m,標準差為3.17 m 的正態分布,墜地撞擊地縱向分布特征如圖10 所示。

圖10 墜地撞擊點縱向分布特征Fig.10 Characteristics of longitudinal distribution of impact points on the ground
按照事故發生危害性評估流程,可計算出初始點的撞擊點分布圖如圖11 所示,計算得出:Sσ= 19.91 m2,S2σ= 79.63 m2,S3σ= 179.17 m2,根據式(11)定義,可計算出風險值為14.81。同理,整個無人機執行任務路徑點均可按照此方式實時計算風險值。事故發生危害性等級標準如表4 所示。對照標準,無人機在初始點的危害等級為輕微危害。

表4 事故發生危害性等級標準Table 4 Accident hazardous level standard

圖11 初始點撞擊點分布圖Fig.11 Distribution of initial point impact points
無人機在運行過程中,通過計算事故危害性數值,結合事故發生危害性等級標準,實時判斷無人機所處的安全水平,整個運行的實時評估結果如圖12 所示。為更精細地顯示無人機整個運行過程中的風險等級,繪制相應的實時風險等級安全譜,如圖13 所示。從安全譜中可以直觀地看出,無人機運行風險在中度危害及以下,只有極個別時刻處于中度危害等級,無嚴重危害等級,整體處于比較安全的水平。

圖12 無人機實時風險評估結果Fig.12 UAV real-time risk assessment results

圖13 無人機實時運行安全譜Fig.13 UAV real-time operation safety spectrum
由此可見,基于提出的風險評估框架能夠實時評估無人機所處的安全水平,為無人機在復雜低空運行提供理論以及與標準。
為驗證本文所提風險評估的優越性,將近年來應用在無人機風險評估領域的安全性分析方法的優缺點進行對比,具體信息如表5 所示。本文所提的風險評估方法能夠無人機墜落的不確定性因素,符合無人機實際飛行情況,更加準確地衡量無人機風險水平。

表5 風險評估方法對比Table 5 Comparison of risk assessment methods
作為UTM 中安全運行能力的體現,風險評估是無人飛行系統能否安全融入國家空域系統的關鍵環節。本文提出考慮無人機墜落過程不確定性的風險評估方法,可將該功能嵌入到UTM 服務中,能夠實時評估無人機的安全水平,用于風險緩解以及突發事件管理。
(1)提出基于無人機墜落動力學模型與Monte-Carlo 方法的事故危害性建模,能夠充分考慮墜落過程中的不確定性,結合事故危害性等級標準,確定無人機運行過程中的風險值。
(2)以108°54'44″N~108°55'13″N,34°14′28″E~34°14'41″E 區域為例,選取大疆經緯 M210-RTK 型無人機作為評估對象,驗證了風險評估方法的可行性。