基于拉格朗日松弛的双线铁路列车运行图优化算法
2016-05-08廖正文苗建瑞孟令云李海鹰
廖正文， 苗建瑞， 孟令云， 李海鹰， 赵 岚
(1. 北京交通大学 轨道交通控制与安全国家重点实验室，北京 100044；2. 北京交通大学 交通运输学院，北京 100044； 3. 西安铁路职业技术学院，陕西 西安 710014)
列车运行图是铁路行车组织的基础，编制和优化运行图是铁路运输组织领域的经典问题。对此，各国专家学者建立目标和约束条件各异的优化模型，并根据模型特点，设计不同的求解算法。但利用数学规划方法对列车运行图进行优化，精确的算法只能求解小规模算例，对于大规模问题必须采取不同的策略将原
问题分解或简化，以求降低问题的规模，如国外CAREY[1]、CAREY and LOCKWOOD[2]、CASTILLO[3]，国内周磊山[4]、彭其渊[5]、史峰[6]等均根据问题的不同特点设计基于问题分解的启发式求解算法。
列车运行图优化，无论是既有线上以旅客列车总旅行时间最短、货物列车开行数量最多，还是客运专线上以列车服务频次高、旅行时间短等为目标，本质均是有限时空资源分配的组合优化问题，尽管当前计算机运算速度和储存空间有较大发展，但是列车运行图的组合规模决定人们依然难以使用传统的方法直接求解大规模问题。因此,目前解决这类问题的思路主要分为两大类：一类为采用模拟人工编图的启发式算法或遗传算法、模拟退火算法等智能优化算法；另一类为将复杂的原问题分解为相对容易求解的子问题，进而获得原问题的较优解。国外有研究采用拉格朗日松弛算法求解列车运行图问题[7-12]，该研究表明拉格朗日松弛算法对于求解列车运行图问题具有一定的适用性，并且能够通过分析优化间隔评价求解结果的质量，是解决列车运行图优化问题的高效方法。
本文针对列车运行图的大规模有限时空资源分配本质问题，结合拉格朗日松弛算法对资源约束问题良好的分解效果及优化程度可评估的特点，以列车对“行车资源”的占用代替传统列车追踪间隔时间表示通过能力约束，考虑实际运营中的列车运行图要素(起停车附加时分、停站时间和各种追踪间隔时间)，构建基于累积流变量的双线铁路列车运行图优化模型；构造拉格朗日松弛算法框架，建立“时间-空间-状态网络”以求解拉格朗日松弛子问题。最后通过实际运营中的列车运行图编制算例验证算法的效率与优化质量。
1 基于累积流变量的列车运行图优化模型
1.1 模型假设和已知条件
为更好地研究影响列车运行图求解算法性能的关键要素，现对列车运行图问题进行适当的抽象和假设。假设列车在双线铁路网上运行，不考虑列车在站内的进路与到发线安排；行车闭塞方式为自动闭塞，追踪间隔时间固定且已知；列车在各区间的纯运行时分、起停附加时分已知且固定；列车的数量、起讫点、径路、停站方案已知；运行图的编制不考虑车底的运用。
模型输入为：列车的等级、不同等级列车的区间纯运行时分、车站起动及停车附加时分、列车追踪间隔时分；列车运行径路，列车始发(接入)时间窗(在始发站或接入站的最早出发时刻至最晚出发时刻之间的时间段)、在各站的最小停站时分以及列车的权重。模型输出为列车在径路上各站的到达和出发时刻。
据文献[13]，若将运行图的时间轴按照最小规划单位离散化，将列车在线路上的运行过程刻画为列车对“行车资源”的占用，利用列车占用“行车资源”的排他性表示列车的追踪间隔，即可采用含有累积流变量的数学等式表达占用“行车资源”的过程，建立数学优化模型，利用拉格朗日松弛方法求解列车运行图问题。因此，本文将利用占用“行车资源”刻画列车的追踪间隔时间，利用累积流变量刻画列车的运行过程。
1.2 行车资源占用模型
对双线铁路，某车站某衔接方向同一时间内最多只能为某一列车办理到达或出发作业，只能依次为各列车单独办理到达或出发进路。从开始办理进路至进路解锁时间段内，其他列车不能办理相应的作业。因此可将车站某衔接方向用于办理到达或出发作业的时间看作待分配的时空资源。本文定义车站上列车占用的“行车资源”为各站列车到达与出发作业在时空上的占用许可，见图1，称为“接车资源”和“发车资源”。列车的到达(或通过)作业需要占用的行车资源自开始办理到达进路时开始，至到达进路解锁时止，见图2(a)。列车出发(或通过)作业需要占用的行车资源自开始办理出发进路时开始，至出发进路解锁时止，见图2(b)。当有列车办理到达(或出发)作业时，接车资源(或发车资源)被占用，否则接车资源(或发车资源)空闲。
传统的列车追踪间隔时间设置目的是为了保证前后行列车占用同一线路资源,在时间上不冲突，追踪间隔时间的值可根据办理列车进路各环节的时间组成计算得到，与本文提出的“行车资源”占用是等价的。因此,可用列车对“行车资源”的占用时间分配刻画列车的追踪间隔，具体的分析为:
列车在车站办理到达或出发进路主要分为3个阶段：
(1) 办理进路、开放信号；
(2) 列车驶过进路；
(3) 列车出清，进路解锁。
各阶段占用时间见图2。
当列车在本站停站时，占用时间根据以下起止时刻确定：
上述占用时空资源要根据列车速度、车站轨道电路结构、接发车进路情况确定。
列车出发追踪间隔时间包括：
(1) 前行列车从车站出发至完全离开车站进入区间的时刻；
(2) 与后行列车办理出发作业时间的组合。
列车到达追踪间隔包括：
(1) 前行列车到达车站时刻起至到达进路解锁的时间；
(2) 办理后行列车到达作业时刻起至该列车到达车站的时间。
列车通过追踪间隔包括：
(1) 前行列车通过车站时刻起至列车完全离开车站进入区间的时刻；
以上分析说明，列车车站追踪间隔时间可用列车占用车站“行车资源”组合等价表示。车站的到达、出发追踪间隔时间一般是列车追踪间隔时间的决定性因素，但对于长大闭塞分区等区间追踪间隔时间是决定性因素的情况，可通过数据预处理，将列车在区间追踪间隔的最大值用该区间前方站的出发间隔时间或后方站的到达间隔时间表示，以统一使用列车的车站行车资源表示追踪间隔。
1.3 基于累积流变量的列车运行过程建模方法
定义路网中车站节点集合为V，区间弧段集合为E，车站节点标号为i,j,k；列车集合F，集合中各列车标号为f。
在列车运行图的时间维度上，采用累积流变量(Cumulative Flow Variables)的建模方法。累积流变量是在时间轴上，将时间按照固定的单位长度进行离散化，并使用一组0-1变量表示某参量在每个离散时刻的状态。
例如图6(b)中，列车在t=5时离开车站i，则t<5时，dfi,t=0，t≥5时，dfi,t=1。在时间轴上表示车站i的累积流变量dfi,t序列时，累积流变量的性质决定该序列最多仅可能在某一时刻发生1次从0到1的突变。在此突变位置t=5之前，即t=1,…,4时，dfi,t=0，dfi,t-1=0，在此突变位置之后，即t=6,7,8,…,T时，dfi,t=1，dfi,t-1=1，该两种情况t×dfi,t-dfi,t-1均为0。只有在突变的位置，即t=5时，dfi,t=1，dfi,t-1=0，则t×dfi,t-dfi,t-1=t。在式(1)求和的各项中，只有在t=5时，该项的值为t，其余项均为0，因此，有等式( 1 )成立。等式( 2 )同理。
( 1 )
( 2 )
( 3 )
( 4 )
列车越行问题本质是列车在区间的运行顺序在某车站发生交换。由于车站接车方向与发车方向的行车资源是分别定义的，二者相互独立。当列车发生越行对于越行站是同一对前后行列车占用接车资源与发车资源的占用顺序发生交换。若车站衔接多方向，可对每个衔接方向的线路各定义一组出发与到达行车资源，当不考虑车站内进路冲突的前提假设下，使模型能够适应多方向情况。
使用累积流变量方法建模可使用线性方程刻画列车对行车资源占用的过程以及行车资源状态随时间变化的过程，即将线性优化与离散系统仿真融合在一起，为轨道交通运输组织优化中引入复杂要素提供手段；其次，用该方法建立的模型，重点在于列车对行车资源的占用关系，而非传统运行图优化模型中的列车与列车运行间隔时间关系，便于对原问题进行分解，实现复杂问题求解的目标。
1.4 基于累积流变量的列车运行图优化模型
对于列车运行图编制的目标，有列车总旅行时间最短、旅客换乘时间最少等。本文研究的重点在于可接受的求解时间内提升列车运行图求解的质量，而导致问题难以获得最优解的原因在于复杂的约束条件，目标函数可以根据不同的编制目标进行改进。因此，选取目标中最具有代表性的加权总旅行时间最短为目标，对于有特定目标的编图问题，目标函数可以在此基础之上扩充。
列车运行图优化模型的目标函数是加权总旅行时间最短为
( 5 )
式中：of、sf为列车f的始发车站和终到站；pf为列车f的权重。由于不同列车等级不同，相同等级的列车也可能存在一定的重要程度差异，此权重体现编图人员对于列车规划优先等级的考虑。
模型约束条件为
(1) 列车到达、出发事件状态约束
dfi,t-1≤dfi,t
∀f∈F∀t∈T:t≠0 ∀i∈V
( 6 )
afi,t-1≤afi,t
∀f∈F∀t∈T:t≠0 ∀i∈V
( 7 )
该组约束保证对于表示列车f在车站上作业的同一组累积流变量随着时间的推移最多只能发生1次变化，且只能从0向1变化，保证列车进入弧段或离开弧段事件从“尚未发生”到“已经发生”的状态转变最多只能发生1次。
(2) 区间进出先后顺序约束
dfi,t≤afj,t
∀f∈F∀t∈T∀(i,j)∈E
( 8 )
该约束保证列车进入区间的时刻不得晚于离开本区间的时刻。
(3) 出发时间窗
dfof,t=0
∀f∈F∀t∈T∧t ( 9 ) dfof,t=1 ∀f∈F∀t∈T∧t>LSTf (10) 式中：ESTf、LSTf分别为列车f最早始发时间和最晚始发时间。 设置出发时间窗约束，可以保证列车在规定的时间范围内出发，尤其是接入到本区段的列车，需要以衔接区段交出列车的时刻作为本区段列车运行图编制的输入条件。同时设置出发时间窗也可以减小问题的解空间，缩小搜索的范围。 (4) 车站最小客运业务停站时间 ∀f∈F∀i∈V (11) 式中：wfi为列车f在i站的最小停站时间。运行图的编制必须满足列车开行方案规定的停靠站点以及停站时间。车站i的出发时刻与到达时刻之差即为列车在该站的停站时间。 (5) 区间运行时分约束 ∀f∈F∀i∈V (12) ∀f∈F∀i∈V (13) FTfi,j+xfi×DAfi+xfj×AAfj= ∀f∈F∀i,j∈E (14) 式(12)～式(14)中：xfi为列车f在i车站是否停车，变量取值为0表示不停车，1为停车；FTf(i,j)为在弧段i,j上的纯运行时分；DAfi、AAf(j)分别为起动附加时分和停车附加时分。 (6) 列车占用行车资源标记约束 ∀f∈F∀t∈T∀i∈V (15) ∀f∈F∀t∈T∀i∈V (16) 式中： (7) 能力约束 每个行车资源被占用的次数是有限的，同一行车资源仅能被一列车占用。 (17) 拉格朗日松弛原理是将问题中的难约束松弛，并将其作为惩罚项添加到目标函数中，使约束条件简化。通过求解拉格朗日松弛问题，可以获得原问题的最优边界，同时通过不断迭代更新拉格朗日乘子，使得松弛问题的解逐渐逼近原问题的最优解[10]。 基于累积变量的列车运行图优化模型中，通常认为制约模型求解效率的约束是能力约束[9]。即式(17)表达的约束。该约束导致列车之间相互作用，从而引发大规模的组合。因此，将其松弛，并将该约束作为惩罚项添加到目标函数中，形成拉格朗日对偶问题Z(D)。 目标函数为 (18) 通过对上述目标函数式进行恒等变换，可得到 (19) 式中： 在拉格朗日对偶问题ZD中，由于松弛约束式(17)，使用式( 6 )～式(16)刻画单个列车运行过程中的时空逻辑，使各列车之间的关联被解除，而目标函数是各列车时空路径的目标函数LRf之和。实质上是列车时空最短路径问题，对于此类问题，通过逐一枚举各列车所有符合约束条件(6)～约束条件(16)可选时空路径，构建时空网络，然后通过最短路径算法快速高效地求得LRf子问题的最优解。 由于约束中包含列车是否停站的决策变量xf(i)，对应的求解子问题的时空网络也可表示列车在车站的作业性质(到达、出发或通过)。因此需要构建改进的列车运行时空网络，将原“时空”网络扩展为“时空状态”网络。时空状态网络有4类节点： (1) 列车的时空逻辑起点和逻辑终点表示列车在时空网络中的生成与消失； (2) 到达节点表示列车的到达作业； (3) 出发节点表示列车的出发作业； (4) 通过节点表示通过列车的到达与出发作业。 在上述的节点上根据列车的作业类型不同构建弧段，正确地表示车站起停附加时分和列车间隔时间。 时空网络中有3类弧段： (1) 起止弧表示列车从逻辑起点进入时空网络，以及列车从时空网络离开前往终止节点的弧段； (2) 运行弧表示列车在两个车站间的运行； (3) 停站弧表示列车在某车站内从到达至出发的状态转变。 通过模型中各种约束建立每列车对应的所有可选时空路径构成1个完整的时空网络，通过求解该时空网络的最短路径可求得问题LRf的最优解，见图8。 在2.1节中，拉格朗日松弛对偶问题ZD在给定拉格朗日乘子时其解是一定的，为使拉格朗日对偶问题的解逐渐逼近原问题的最优解，要通过迭代更新拉格朗日乘子的值。迭代过程中，采用次梯度方法对拉格朗日乘子进行更新，为 (20) 式(20)表示新的拉格朗日乘子与当前的拉格朗日乘子及本次迭代求解得到的时空资源占用次数有关，在当前拉格朗日乘子的基础上，根据本次求解中列车对该资源的占用次数，对拉格朗日乘子进行调整。 在拉格朗日松弛问题中，拉格朗日乘子的现实含义是“时空资源费用”，即列车占用“时空资源”的代价，更新拉格朗日乘子的次梯度法则是根据本次迭代中列车对资源的竞争结果情况对“时空资源费用”进行动态更新。将列车对时空资源的请求看作“需求”，将线路能力看作“供给”，将拉格朗日乘子的值看作“价格”。各次列车首先根据时间效益最大的路径，独立地向时空资源发出占用请求，各时空资源根据所有列车的占用请求与自身能力的匹配情况，使用次梯度法更新拉格朗日乘子值，即根据“供需关系”对各资源“定价”，从而影响下次迭代中各列车在时空网络中的路径选择。通过“资源定价”与“路径选择”过程的不断迭代，从而解决运行图优化的核心问题——有限时空资源的优化分配。 由于使用拉格朗日松弛方法将原问题分解后，得到的解可能是不可行解。因此需要提供高效获得原问题可行解的方法。现使用启发式算法，即把拉格朗日松弛求解得到的结果作为启发信息，按照一定的顺序逐个列车进行规划，从而求得可行解。 该启发式方法在本质上是滚动式求解方法，关键是如何确定进入规划列车的顺序。采用按照列车的LRPf由大到小的顺序对列车进行滚动优化求解，LRPf的计算方法为 LRPf= (21) LRPf的值越大说明该列车冲突越严重，冲突疏解的难度较大，可调整余地较小；该值越小说明该列车冲突越少，可调整余地越大，故按照由难到易的顺序逐个对列车的运行线进行铺画。 铺画过程中，使用与求解拉格朗日松弛问题同样的最短路径方法进行求解，并同时对占用进行标记。算例验证中，如果该启发式方法不能获得可行解，即有列车不能被规划，则调整这些规划失败列车的优先级，优先铺画前次规划失败的列车。 将京广高速铁路武汉—广州南区段(共18座车站)作为算例对模型与算法的可行性和求解效率进行分析。使用ILOG CPLEX V12.3软件求解原问题与线性松弛问题，拉格朗日松弛算法框架在Visual Studio 2012平台上使用C#语言实现，算法运行环境为1台CPU为Intel Xeon E5520 2.27 GHz×2，8 GB内存的台式计算机。 首先采用武汉—长沙南区段7座车站对算法的可行性和性能进行初步验证。算例设置求解总时间长度为180 min，时间离散化后的铺画时间单位为1 min，列车出发时间窗为全时间域(0～180 min)，列车权重均设为相等，测试随着列车数量增加求解质量与求解时间的变化。实验首先使用CPLEX分别计算模型的最优解与线性松弛解(下界)。针对相同实验数据集，再次利用拉格朗日松弛算法框架进行200次迭代计算。求解结果表明，随着列车数量增大，CPLEX求最优解变得越来越困难，当列车数增加到8列时，在设定时间内(7 200 s)不能完成求解过程，但将原问题进行线性松弛后，CPLEX能够较快地求解出问题的下界。 表1 算法求解性能分析 与CPLEX求解出的线性松弛生成的下界相比，求解8列车时，拉格朗日松弛算法得到的下界质量有11.13%的提升。 进一步分析CPLEX与拉格朗日松弛算法在求解问题下界方面的效率，当列车数量较少时，CPLEX求解速度比使用拉格朗日松弛算法快，但是当列车数量增多(大于32列)时，拉格朗日松弛算法具有明显的效率优势。 使用拉格朗日松弛算法时，上下界间隔将随着列车的数量增加而增加，见图9。 以武汉—广州南全线下行方向为背景测试模型和算法的性能。设定全天(1 440 min)下行方向开行110列车，时间离散化后的铺画时间单位为1 min，为每列车设置30 min的始发站出发时间窗，并根据列车的等级、停站数量、旅行时间设置权重。利用拉格朗日松弛算法进行200次迭代，总耗时10 448 s，收敛的上下界间隔为0.19%，见图10。 根据我国双线铁路列车运行图编制特点，考虑列车起停附加时分和列车追踪间隔等关键要素，构建基于累积流变量的列车运行图0-1整数规划模型。 针对运行图优化问题求解困难的问题，结合模型特点提出拉格朗日松弛算法框架，对问题按照列车分解，针对子问题特点，提出改进的列车运行时空状态网络并设计子问题的最短路径求解方法。 算例结果表明，拉格朗日松弛算法能够高效地求解大规模列车运行图问题，并能够通过比较优化上下界的间隔确定运行图的优化程度。同时说明，拉格朗日松弛方法所求得下界质量比线性松弛更优，使用启发式方法求解上界速度比较快。当列车数量增加时，运行线之间的冲突将更加频繁。因此，可能导致拉格朗日松弛算法的上下界间隔随着列车的数量增加而逐渐变大。 参考文献： [1] CAREY M. A Model and Strategy for Train Pathing with Choice of Lines，Platforms，and Routes[J]. Transportation Research Part B：Methodological，1994，28(5)：333-353. [2] CAREY M，LOCKWOOD D. A Model，Algorithms and Strategy for Train Pathing[J]. Journal of the Operational Research Society，1995，46(8)：988-1 005. [3] CASTILLO E，GALLEGO I，URENA J M，et al. Timetabling Optimization of A Mixed Double and Single-tracked Railway Network[J]. Applied Mathematical Modelling，2011，35(2)：859-878. [4] 周磊山，胡思继，马建军，等. 计算机编制网状线路列车运行图方法研究[J]. 铁道学报，1998，20(5)：16-22. ZHOU Leishan，HU Siji，MA Jianjun，et al. Network Hierarchy Parallel Algorithm of Automatic Train Scheduling[J]. Journal of the China Railway Society，1998，20(5)：16-22. [5] 彭其渊，朱松年，王培. 网络列车运行图的数学模型及算法研究[J]. 铁道学报，2001，23(1)：1-8. PENG Qiyuan，ZHU Songnian，WANG Pei. Study on A General Optimization Model and Its Solution for Railway Network Train Diagram[J]. Journal of the China Railway Society，2001，23(1)：1-8. [6] 史峰，黎新华，秦进，等. 单线列车运行图铺划的时间循环迭代优化方法[J]. 铁道学报，2005，27(1)：1-5. SHI Feng， LI Xinhua，Qin Jin，et al. A Timing-cycle Iterative Optimizing Method for Drawing Single-track Railway Train Diagrams[J]. Journal of the China Railway Society，2005，27(1)：1-5. [7] CAPRARA A，FISCHETTI M，TOTH P. Modeling and Solving the Train Timetabling Problem[J]. Operations Research，2002，50(5)：851-861. [8] CAPRARA A，MONACI M，TOTH P. A Lagrangian Heuristic Algorithm for A Real-world Train Timetabling Problem[J]. Discrete Applied Mathematics，2006，154(5)：738-753. [9] MENG Lingyun，ZHOU Xuesong. Simultaneous Train Rerouting and Rescheduling on An N-track Network：A Model Reformulation with Network-based Cumulative Flow Variables[J]. Transportation Research Part B：Methodological，2014，67：208-234. [10] FISHER M L. An Applications Oriented Guide to Lagrangian Relaxation[J]. Interfaces，1985，15(2)： 10-21. [11] GOH C J，MEES A I. Optimal Control on A Graph with Application to Train Scheduling Problems[J]. Mathematical and Computer Modelling，1991，15(2)：49-58. [12] BRANNLUND U，LINDBERG P O，NOU A，et al. Railway Timetabling Using Lagrangian Relaxation[J]. Transportation Science，1998，32(4)：358-369. [13] HARROD S. A Tutorial on Fundamental Model Structures for Railway Timetable Optimization[J]. Surveys in Operations Research and Management Science，2012，17(2)：85-96.
2 基于拉格朗日松弛的模型求解方法
2.1 拉格朗日松弛问题
2.2 基于时空状态网络最短路径的子问题求解方法
2.3 拉格朗日乘子更新方法
2.4 获得可行解的启发式算法
3 算例分析
4 结论
