CN103226900B - 一种基于加权图模型的空域扇区划分方法 - Google Patents

Abstract

本发明属于空中交通管理领域，提出了一种基于加权图模型的空域扇区设计方法。本发明在建立一种能够准确表示航路和空中交通量的无向图模型后，将图顶点作为Voronoi图基点把空域离散化，然后依据每个Voronoi图单元的工作负荷和航路上的交通量构建加权图模型。进而，利用融合了一般加权图切算法、负荷平衡算法和启发式算法的图划分方法，将加权图模型划分为多个子图，再由各个子图所包含的顶点映射到对应的Voronoi图单元组合形成扇区。由本方法设计的扇区不仅满足负荷平衡和协作负荷最小化约束，而且满足扇区最小距离约束、凸性约束和连通性约束。

Description

一种基于加权图模型的空域扇区划分方法

技术领域

本发明属于空中交通管理技术领域，涉及一种空中交通管理中的空域扇区设计方案，尤其是一种能够满足扇区负荷均衡和交接负荷最小化以及能够保证扇区凸性、最小距离和连通性等几何约束的扇区划分方法。

背景技术

目前，世界上绝大多数国家的空域组织形式是固定航路和固定扇区，即飞机沿着固定航路飞行,扇区边界和扇区数目是固定不变的。这是一种结构化的静态空域。这种静态空域在恶劣天气和空中交通拥挤的情况下容易造成航班延误。为了解决上述交通难题，研究者提出了满足空域使用者需要的动态空域配置概念。在动态配置下，飞机可以按照适合自己的航线自由飞行，从而使得航路结构从静态转向动态。这种动态交通模式对现有的空中交通管理特别是现有的扇区结构是一种挑战，因为目前扇区结构是静态的，是在一定历史条件下依靠经验形成的，已经无法满足动态空域配置下交通管理的需要，因此需要按照动态变化的空中交通量来设计空域扇区。

近年来，有关学者针对空域扇区设计做了一些有益的研究和探索。从现有文献可知绝大多数扇区划分方法是首先建立空域模型，然后按照管制员负荷均衡原则利用优化算法将空域模型划分成扇区。所用空域模型分为：细胞元模型、航路模型、Voronoi图模型和图模型。利用细胞元模型，给定空域被划分成一系列正六边形网格，然后利用优化算法将这些正六边形网格组合形成扇区。例如，Yousefi于2004年发表的《Temporal and Spatial Distribution ofAirspace Complexity for New Methodologies in Airspace Design》，Klein于2005年发表的《AnEfficient Method for Airspace Analysis and Partitioning based on Equalized Traffic Mass》，Drew于2008年发表的《Analysis of an Optimal Sector Design Method》和Tien于2009年发表的《Optimizing Airspace Sectors for Varying Demand Patterns using Multi-Controller Staffing》分别利用厂址选择算法，区域生长算法数规划算法将细胞元模型组合形成扇区。然而在上述文献中所设计的扇区边界是锯齿状，扇区在形状上并不尽人意。利用航路模型，Briton于2008年发表的《Airspace Partitioning Using Flight Clustering and computational Geometry》将航路点利用k-means算法聚类成扇区，而Basu于2008年发表的《Geometric Algorithms for OptimalAirspace Design and Air Traffic Controller Workload Balancing》则发展了一种几何算法来解决扇区设计问题。使用Voronoi图模型，Delahaye于1998年发表的《Airspace Sectoring byEvolutionary Computation》建立了初始扇区，然后利用遗传算法通过移动扇区边界以达到负荷均衡。更进一步地，Xue于2009年发表的《Airspace Sector Redesign Based on VoronoiDiagrams》在Delahaye扇区划分方法基础上，应用迭代加深算法来改进遗传算法在实现上的不足。然而，基于上述三种模型的扇区划分方法没有考虑到空域结构因素，因而所设计的扇区并不能完全保证扇区的凸性、最小距离和连通性等几何约束。而图模型通过则充分考虑了空域结构信息，将机场、航路点和航路交叉点等关键点用图的顶点表示，航路用边来表示，建立空域的图模型，进而在图模型基础上方便地采取措施进行模型优化，从而获得扇区。例如，Trandac于2003年发表的《Optimized Sectorization of Airspace with Constraints》在图模型基础上利用约束算法优化形成扇区。Li于2010年发表的《Graph-Based Algorithm for DynamicAirspace Configuration》则首先建立一个边被赋予了交接负荷的图模型，然后利用谱聚类对图模型进行划分，但这样并不能保证子图的工作负荷均衡，因而又通过合并相邻子图的方式以获得负荷均衡且满足连通性约束的子图，最后通过确定边界形成扇区。而Martine于2007年发表的《A Weighted-graph Approach for Dynamic Airspace Configuration》将顶点和边都没有赋予工作负荷的图模型，利用谱二分法将按照所包含的航路点和机场数目均衡以及连接两部分的航路条数最小化原则将网格模型迭代地进行二分成许多子图。然而，由于空域的交通量分布并不完全一致，所以扇区间航路最小化并不能扇区间交接负荷的最小化，航路点和机场数目的均衡并不能保证工作负荷的均衡，所以又将空域离散成众多正四方形网格，通过子图间网格调整获得负荷均衡的扇区。即便如此，也不能保证扇区间交接负荷的最小化。而且，利用正方形网格也不能完全保证扇区最小距离约束和凸性约束。

综上所述，现有的方法绝大多数方法并不能完全保证扇区的几何约束，而且在扇区设计时未能兼顾扇区的负荷均衡和交接负荷最小化。

发明内容

为解决现有技术中存在的几何约束不能得到完全保证及未能兼顾负荷均衡和交接负荷最小化的问题，本发明提供一种基于加权图模型的扇区划分方法，由本方法设计的扇区不仅能够满足扇区之间负荷平衡和交接负荷最小化，而且能够满足扇区的凸性、最小距离和连通性约束等几何约束。

本发明采用如下的技术方案。

一种基于加权图模型的空域扇区划分方法，其特征在于包括以下步骤：

步骤一，建立空域加权图模型，方法如下：

(1)根据给定的空域结构信息，构建无向图模型G＝G(V,E)，其中顶点集V＝{1,2,...,n}，代表机场、航路点和航路交叉点等关键点，边集合E＝{e_ij＝(i,j):i,j∈V}，e_ij表示连接顶点i和j的航路段。

(2)在无向图模型基础上，以顶点集合V中的顶点作为基点建立Voronoi图D。利用D将空域划分为n个Voronoi图单元D_i(i=1,2，...,n)，其中D_i表示第i个顶点所对应的Voronoi图单元。Voronoi图单元的边界在合并单元成为扇区时将部分地成为扇区边界。这种通过Voronoi图离散空域模型形成单元的方式可以保证将单元组合形成扇区时满足扇区的凸性约束。对空域离散化后，有一些基点离单元边界很近，如果所设计的扇区边界恰好就落在这条边界上，那么扇区就不符合最小距离约束，所以将这条边界删除，将此条边界所连接的两个单元合并。合并后的单元数目为r(r≤n)，通过合并措施使得这些单元在组合成扇区时能够保证最小距离约束。

(3)在空域离散化的基础上，根据交通数据计算一段时间内D_i(i=1,2，...,r)的工作负荷w_i和每条航路e_ij上的交接负荷。

(4)为了便于分析，本发明将D_i抽象成一个辅助的加权图模型，其中D_i抽象成一个顶点，顶点权重就是D_i的工作负荷w_i；若D_i和D_j之间有航路连接，则所有航路就抽象成一条边，边的权重就是连接D_i和D_j所有航路上交接负荷之和，记为w_ij。将所构造的加权图模型记为G_w，G_w=(w,W)，其中w为顶点权重向量，W为边的权重矩阵表示：

w＝[w₁,w₂,...,w_r]^T     (1)

W＝[w_ij]_r×r，w_ij=w_ji  (2)

至此，扇区划分问题就转换为一个图划分问题，所构造的加权图模型顶点和Voronoi图模型单元之间存在着一一对应关系，这种对应关系可以方便地将对加权图操作时所得到的性质映射到Voronoi图中。

步骤二，确定扇区数目，方法如下：

根据给定空域的总负荷A_c和每个扇区的最大负荷S_c确定扇区数k，也就是确定划分加权图子图数目，扇区数k按下式计算：

其中，是向上取整运算。

步骤三，划分加权图模型：本发明发展了一种由一般加权图切算法、负荷平衡算法和由K-L算法增益概念所得到的启发式算法融合成的图划分方法，对空域加权图模型进行划分。图划分方法分为两步：一是将图划分为k个子图，二是利用负荷平衡算法和启发式算法在子图间移动顶点以实现负荷平衡。具体方法如下：

(1)利用一般加权图切算法，按照下式将加权图模型划分为k个子图(i=1,2,…,k)：

$\underset{G_w^1,...,G_w^k}{\min }\underset{i=1}{\overset{k}{\mathrm{\Σ}}}\frac{\mathrm{cut}(G_w^i,G_w/G_w^i)}{w\left(G_w^i\right)}$

$w\left(G_w^i\right)=\underset{v_c\∈G_w^i}{\mathrm{\Σ}}{w}_c$

$\mathrm{cut}(G_w^i,G_w/G_w^i)=\underset{v_c\∈G_w^i,v_d\∉G_w^i}{\mathrm{\Σ}}{w}_{\mathrm{cd}}$

其中，V_c、V_d为顶点，w_c为顶点V_c的权重，w_cd为加权图边的权重；

1）构造如(1)和(2)式所示的顶点权重向量w和协作负荷矩阵W；

2）求W每列之和，将该和作为对角元素构造对角矩阵D，按下式计算拉普拉斯矩阵L：

L＝D-W

3）计算

$s\left(i\right)=\frac{1}{\sqrt{w\left(i\right)}}(i=\mathrm{1,2},...,r)$

由s(i)构成对角矩阵：

S＝diag{s(1),s(2),...,s(r)}

4)按下式构造矩阵C：

C＝S×L×S

5）计算C的特征值和特征向量，将特征值按照从小到大的顺序排列，取前k个特征值对应的特征向量，利用k-means聚类算法将加权图分为k类。

因为子图顶点权重之和为：

$w\left(G_w^i\right)=\underset{v_d\∈G_w^i}{\mathrm{\Σ}}{w}_d$

其中，w_d为顶点V_d的权重；

连接子图的边权重之和为：

$\underset{i=1}{\overset{k}{\mathrm{\Σ}}}\mathrm{cut}(G_w^i,G_w/G_w^i)=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}\underset{v_d\∉G_w^i}{\underset{v_c\∈G_w^i,}{\mathrm{\Σ}}}{w}_{\mathrm{cd}}$

所以，切算法兼顾了子图顶点权重之和的均衡与连接子图的边权重之和的最小化，也就是说，在连接子图的边权重之和最小化条件下的子图顶点权重之和并不一定完全相等，因而需要采用负荷平衡算法和启发式算法进行子图间的顶点移动以达到子图间的顶点权重之和相等，其中负荷平衡算法用来计算移动顶点权重多少，启发式算法用来移动顶点。

(2)采用负荷平衡算法计算需要移动的顶点权重数目，算法实现如下：

1)按下式计算子图平均权重和向量b：

$\stackrel{\‾}{w}=\frac{1}{k}\underset{i=1}{\overset{k}{\mathrm{\Σ}}}w\left(G_w^i\right)$

$b={(w\left(G_w^1\right)-\stackrel{\‾}{w},w\left(G_w^2\right)-\stackrel{\‾}{w},...,w\left(G_w^k\right)-\stackrel{\‾}{w})}^T$

2)构建子图的关联矩阵Z，按下式计算拉普拉斯矩阵L：

L＝ZZ^T

3）构建方程式Ld＝b，解出d，按下式求子图间需要移动权重x_ij：

x_ij＝Zd^T

(3)采用启发式算法移动顶点，获得一系列满足顶点权重均衡条件下连接边权重之和最小化的子图(i=1,2,…,k)，方法如下：

首先，确定子图和之间顶点的移动方向。假定顶点将从移动到中，确定中与相邻的顶点集B_ij，计算B_ij中所有顶点的权重之和a_ij和每个顶点的增益g_f：

$g_f=\underset{v_d\∈G_w^i,v_f\∈G_w^j}{\mathrm{\Σ}}{w}_{\mathrm{df}}-\underset{v_d,v_e\∈G_w^i}{\mathrm{\Σ}}{w}_{\mathrm{de}}$

其中，v_d、v_e、v_f为顶点，w_df和w_de为加权图边的权重。

然后，将顶点按照增益自大到小降序排列，将增益最大的顶点先移动到如果顶点权重小于x_ij，则移动增益第二大的顶点，重复这个步骤移动顶点。如果a_ij小于x_ij，在移动完B_ij中所有顶点后，在中会有新的顶点毗邻重复上述计算增益和移动顶点的步骤，直至满足x_ij要求。这样，就得到了一系列满足在顶点权重均衡的条件下连接边权重之和最小化的子图(i=1,2,…,k)，即：

$\underset{G_w^{1^{\′}},...,G_w^{k^{\′}}}{\min }\underset{c=1}{\overset{k}{\mathrm{\Σ}}}\frac{\mathrm{cut}(G_w^{i^{\′}},G_w/G_w^{i^{\′}})}{w\left(G_w^{i^{\′}}\right)}$

$w\left(G_w^{i^{\′}}\right)=w\left(G_w^{j^{\′}}\right)$

其中， $w\left(G_w^{i^{\′}}\right)=\underset{v_c\∈G_w^{i^{\′}}}{\mathrm{\Σ}}{w}_c,$ $\mathrm{cut}(G_w^{i^{\′}},G_w/G_w^{i^{\′}})=\underset{v_d\∉G_w^{i^{\′}}}{\underset{v_c\∈G_w^{i^{\′}},}{\mathrm{\Σ}}}{w}_{\mathrm{cd}},$ v_c、v_d为顶点，w_cd为加权图边的权重。

步骤四，将每个子图所包含的顶点对应的Voronoi图单元组合形成扇区。根据加权图模型顶点和Voronoi图模型单元之间的对应关系，所划分的扇区能够满足负荷均衡和交接负荷最小化。

本发明的有益效果是，所设计的扇区不仅能够均衡地分配工作负荷，而且穿越扇区间的交接负荷最小化，这样可以尽可能地减小管制员工作强度。同时，保证了扇区凸性、最小距离和连通性约束，使得管制员有足够的时间去解决飞机飞行期间可能发生的各种冲突，提高了飞机飞行的安全性。

附图说明

图1为基于加权图模型的空域扇区划分方法流程图；

图2为构建加权图模型的流程图；

图3为空域无向图模型；

图4为以图3的顶点为基点建立的Voronoi图；

图5为对关键点距单元边界过近情况进行预处理的Voronoi图；

图6为构建的加权图模型；

图7为划分加权图模型的流程图；

图8为利用一般加权图切算法对图模型进行划分的结果示意图；

图9为顶点移动后的结果示意图；

图10为子图顶点映射到Voronoi图单元形成的扇区划分结果图；

图11为北京管制区现在用扇区划分图；

图12为北京管制区现有扇区与本发明扇区划分方案划分结果对比曲线。

具体实施方式

下面结合附图和实施例对本发明作进一步说明。本实施例针对北京管制区利用真实的交通数据进行空域扇区划分。

图1是本发明所述的方法流程图，具体包括以下步骤：

步骤一，构建加权图模型，方法如图2所示，包括以下内容：

(1)根据给定空域结构信息给出空域的无向图模型，如图3所示，顶点表示机场、航路交叉点和航路点等关键点，顶点间的线段代表航路。

(2)以图3的顶点作为基点，构建如图4所示的Voronoi图，其中细线表示航路，粗线表示Voronoi图边界。因为有些单元的基点距离单元边界过近，所以对单元进行预处理，处理后的结果如图5所示。

(3)计算各个单元在13:00-15:00内的工作负荷和航路上的交接负荷。

(4)将Voronoi图单元抽象成一个加权图模型G_w，构建的加权图模型如图6所示。

步骤二，确定扇区数目，方法如下：

根据给定空域的总负荷A_c（18）和每个扇区的最大负荷S_c（8）确定扇区数k，即划分加权图子图数目：

所以，需要将空域划分为3个扇区。

步骤三，进行加权图划分，方法如图7所示，包括以下内容：

(1)利用一般加权图切算法对加权图进行划分，划分结果如图8所示，虚线表示子图分界线。

(2)采用负荷平衡算法计算需要移动的顶点权重数目。根据切算法划分结果，可知 故需要从子图3中移动负荷x₂₃给子图4。由平衡算法可知x₂₃＝6。

(3)利用启发式算法进行顶点移动，移动后的子图划分如图9所示。

步骤四，将三个子图中分别包含的顶点映射到Voronoi图中所对应的单元，就形成了3个扇区，如图10所示。

应用本发明所划分的扇区如附图10所示，现有扇区划分如附图11所示，二者的对比结果如图12所示。从图12可以看出：

(1)本发明所划分的扇区为3个，扇区负荷分别为5、6、7。而现有扇区为6个，负荷为4、2、0、3、4、5。本发明所设计的扇区不仅少于现有扇区，而且工作负荷更加均衡。

(2)本发明所划分的扇区的交接负荷为13，而现有扇区的为60个。这说明本发明所划分的扇区交接负荷要远小于现有扇区的交接负荷，从而管制员在管理同样数目的飞机时有更小的工作压力。

(3)本发明所划分的扇区边界满足最小距离、连通性和凸性约束，而现有扇区的边界有的恰好在航路交叉点上，不符合最小距离约束。

Claims (2)

1.一种基于加权图模型的空域扇区划分方法，其特征在于包括以下步骤：

步骤一，建立空域加权图模型，方法如下：

(1)根据给定的空域结构信息，构建无向图模型G＝G(V,E)，其中，顶点集V＝{1,2,...,n}代表机场、航路点和航路交叉点等关键点，边集合E＝{e_ij＝(i,j):i,j∈V}，e_ij表示连接顶点i和j的航路段；

(2)在无向图模型基础上，以顶点集合V中的顶点作为基点建立Voronoi图D；利用D将空域划分为n个Voronoi图单元D_i，i＝1,2，...,n，其中D_i表示第i个顶点所对应的Voronoi图单元；Voronoi图单元的边界在合并单元成为扇区时将部分地成为扇区边界；这种通过Voronoi图离散空域模型形成单元的方式可以保证将单元组合形成扇区时满足扇区的凸性约束；对空域离散化后，有一些基点离单元边界很近，如果所划分的扇区边界恰好落在这条边界上，那么扇区不符合最小距离约束，将这条边界删除，将此条边界所连接的两个单元合并；合并后的单元数目为r，r≤n，通过合并措施使得这些单元在组合成扇区时能够保证最小距离约束；

(3)在空域离散化的基础上，根据交通数据计算一段时间内D_i的工作负荷w_i和每条航路e_ij上的交接负荷，i＝1,2，...,r；

(4)将D_i抽象成一个辅助的加权图模型，其中，D_i抽象成一个顶点，顶点权重就是D_i的工作负荷w_i；若D_i和D_j之间有航路连接，则所有航路就抽象成一条边，边的权重就是连接D_i和D_j所有航路上交接负荷之和，记为w_ij；将所构造的加权图模型记为G_w，G_w＝(w,W)，其中w为顶点权重向量，W为边的权重矩阵表示：

w＝[w₁,w₂,...,w_r]^T         (1)

W＝[w_ij]_r×r，w_ij＝w_ji         (2)

步骤二，确定扇区数目，方法如下：

根据给定空域的总负荷A_c和每个扇区的最大负荷S_c确定扇区数k，也就是确定划分加权图子图数目，扇区数k按下式计算：

其中，是向上取整运算；

步骤三，采用由一般加权图切算法、负荷平衡算法和启发式算法融合成的图划分方法划分加权图模型；

步骤四，将每个子图所包含的顶点对应的Voronoi图单元组合形成扇区。

2.根据权利要求1所述的一种基于加权图模型的空域扇区划分方法，其特征在于，步骤三划分加权图模型的方法包括以下步骤：

(1)利用一般加权图切算法，按照下式将加权图模型划分为k个子图i＝1,2,…,k：

$\begin{array}{cc}\underset{G_w^1,...,G_w^k}{\min }& \underset{c=1}{\overset{k}{\mathrm{\Σ}}}\frac{\mathrm{cut}(G_w^i,G_w/G_w^i)}{w\left(G_w^i\right)}\end{array}$

$w\left(G_w^i\right)=\underset{v_c\∈G_w^i}{\mathrm{\Σ}}{w}_c$

$\mathrm{cut}(G_w^i,G_w/G_w^i)=\underset{v_c\∈G_w^i,v_d\∉G_w^i}{\mathrm{\Σ}}{w}_{\mathrm{cd}}$

其中，V_c、V_d为顶点，w_c为顶点V_c的权重，w_cd为加权图边的权重；

1)构造如(1)式和(2)式所示的顶点权重向量w和协作负荷矩阵W；

2)求W每列之和，将该和作为对角元素构造对角矩阵D，按下式计算拉普拉斯矩阵L：

L＝D-W

3)计算

$s\left(i\right)=\frac{1}{\sqrt{w\left(i\right)}},i=\mathrm{1,2},...,r$

由s(i)构成对角矩阵：

S＝diag{s(1),s(2),...,s(r)}

4)按下式构造矩阵C：

C＝S×L×S

5)计算C的特征值和特征向量，将特征值按照从小到大的顺序排列，取前k个特征值对应的特征向量，利用k-means聚类算法将加权图分为k类；

(2)采用负荷平衡算法计算需要移动的顶点权重数目，算法实现如下：

1)按下式计算子图平均权重和向量b：

$\stackrel{\‾}{w}=\frac{1}{k}\underset{i=1}{\overset{k}{\mathrm{\Σ}}}w\left(G_w^i\right)$

$b{(w\left(G_w^1\right)-\stackrel{\‾}{w},w\left(G_w^2\right)-\stackrel{\‾}{w},...,w\left(G_w^k\right)-\stackrel{\‾}{w})}^T$

2)构建子图的关联矩阵Z，按下式计算拉普拉斯矩阵L：

L＝ZZ^T

3)构建方程式Ld＝b，解出d，按下式求子图间需要移动权重x_ij：

x_ij＝Zd^T

(3)采用启发式算法移动顶点，获得一系列满足顶点权重均衡条件下连接边权重之和最小化的子图i＝1,2,…,k，方法如下：

首先，确定子图和之间顶点的移动方向；假定顶点将从移动到中，确定中与相邻的顶点集B_ij，计算B_ij中所有顶点的权重之和a_ij和每个顶点的增益g_f：

$g_f=\underset{v_d\∈G_w^i,v_f\∈G_w^j}{\mathrm{\Σ}}{w}_{\mathrm{df}}-\underset{v_d,v_e\∈G_w^i}{\mathrm{\Σ}}{w}_{\mathrm{de}}$

其中，v_d、v_e、v_f为顶点，w_df和w_de为加权图边的权重；

然后，将顶点按照增益自大到小降序排列，将增益最大的顶点先移动到如果顶点权重小于x_ij，则移动增益第二大的顶点，重复这个步骤移动顶点；如果a_ij小于x_ij，在移动完B_ij中所有顶点后，在中会有新的顶点毗邻重复上述计算增益和移动顶点的步骤，直至满足x_ij要求；这样，就得到了一系列满足在顶点权重均衡的条件下连接边权重之和最小化的子图i＝1,2,…,k，即：

$\begin{array}{cc}\underset{G_w^1,...,G_w^k}{\min }& \underset{c=1}{\overset{k}{\mathrm{\Σ}}}\frac{\mathrm{cut}(G_w^{i^{\′}},G_w/G_w^{i^{\′}})}{w\left(G_w^{i^{\′}}\right)}\end{array}$

$w\left(G_w^{i^{\′}}\right)=w\left(G_w^{j^{\′}}\right)$

其中， $w\left(G_w^{i^{\′}}\right)=\underset{v_c\∈G_w^{i^{\′}}}{\mathrm{\Σ}}{w}_c,\mathrm{cut}(G_w^{i^{\′}},G_w/G_w^{i^{\′}})=\underset{\underset{v_d\∉G_w^{i^{\′}}}{v_c\∈G_w^{i^{\′}}},}{\mathrm{\Σ}}{w}_{\mathrm{cd}},$ v_c、v_d为顶点，w_cd为加权图边的权重。