基于Pareto多目标优化的高速铁路列车节能驾驶曲线集生成

张惠茹，贾利民，王 莉

(1.北京交通大学 轨道交通控制与安全国家重点实验室，北京 100044； 2.北京交通大学 交通运输学院，北京 100044；3.北京市城市交通信息智能感知与服务工程技术研究中心，北京 100044; 4.城市轨道交通系统安全与运维保障国家工程实验室，广东 广州 510330)

随着铁路运营里程的不断增长和运行速度的不断提升，铁路领域的可持续发展受到越来越多的关注，而节能驾驶是目前最直接的有效举措。实际的列车操纵是一个能耗与时间相博弈的过程，即多目标优化问题。目前国内外相关学者在列车节能驾驶模式和多目标优化等方面做了积极的探索。

文献[1]最早研究了列车节能驾驶策略，并利用Pontryagin最大原则求解限速条件下的最优控制问题；文献[2]进一步推导出地铁系统的最优控制序列，不含匀速阶段；文献[3]利用Hamiltonian函数推导出匀速运行以及惰行点是列车节能的关键。以上文献所述为拥有严谨理论推导的第一类方法，但需对列车、线路等条件进行简化，而如果考虑实际运行环境引入复杂的非线性方程或约束，则问题的简化会很困难。文献[4]将基本的列车模型改写为以动能和时间为状态量的连续空间模型，并运用分段仿射近似为混合整数规划求解。这种方法为第二类算法，即将原始问题转化为数学优化问题然后直接求解，但是该类方法建立在大量的假设基础上，且对模型精度要求高。文献[5]提出一种通过仿真分配列车富余时间以达到节能控制的方法；文献[6]设计了一种基于启发式遗传算法的列车牵引运行算法。这两种方法为结合仿真或进化算法的第三类方法，得益于计算机技术的发展，列车、线路、牵引模式等复杂的现实条件可以很容易地加入到模型中，从而得到更加符合实际情况的有效的驾驶模式。

对于所得的驾驶策略，除了赋予目标不同权重得到一个最优解[5]，也可将列车节能驾驶抽象为多目标优化问题获得一组解，而后者为铁路管理者提供了可行的备选方案。基于地铁列车自动驾驶系统提供的加速度、减速度、匀速度以及惰行速度的范围，文献[7]提出了一种多目标粒子群优化算法，并基于仿真得到了Pareto前沿。文献[8]设计一种多目标混合整数精英遗传算法，在离散化轨道线路的基础上，获得了地铁系统能耗与时间的Pareto曲线。文献[9]针对高速铁路系统的节能驾驶问题，设计了一种结合仿真的多目标求解算法，但利用平均等效坡度对线路条件进行了较大程度的简化

综上，既有研究多以地铁为研究对象，操纵策略较为简单，且很少考虑实际的复杂线路条件；部分涉及高速铁路的研究中，也往往将复杂的模型参数进行简化处理，以便快速求得单一的列车运行操纵策略。本文在考虑高速铁路实际线路条件基础上，重点面向列车运行环境复杂多样、易受不同程度的扰动，从而偏离基本计划的现实问题，建立以能耗、时间为目标，以列车受力、区间限速等为约束的多目标列车操纵优化模型，同时设计复杂线路条件简化方法以及改进NSGA-Ⅱ算法，并结合变步长仿真，实现多目标节能驾驶曲线集的生成。该研究结果可针对不同运行条件提供丰富的列车节能操纵策略。

1 问题描述

高速铁路相邻站间线路距离长，且线型条件复杂，包括频繁变化的坡度、曲线、线路限速等。在进行牵引计算时，可以根据线路条件划分计算区段，即一个计算区段内线路条件相同。本文采用最大加速度加速、无制动调速，以及停车前最大减速度制动三种节能策略[10]。其中，在加速阶段使用最大加速度、制动阶段使用最大减速度已经被证明是有效的节能策略[1]。考虑到制动过程将造成能耗浪费，节能驾驶策略应尽量避免非必要的制动，因此采用无制动调速，即在非停车计算区段，列车采用牵引或惰行工况调节速度的策略。

如图1所示，根据坡度、曲线和限速条件，线路分为4个计算区段(Seg1—Seg4)，由相邻计算区段的速度连续关系，即可得到整个区间的速度-距离曲线，而列车节能的关键在于不同计算区段内驾驶策略的选择。

图1 一种节能驾驶的速度-距离曲线示意

2 模型建立

基于列车运行特性及本文所采用的节能策略，在满足受力、速度、距离等约束条件下，建立多目标优化模型。

2.1 假设与变量

为方便描述问题，作以下假设：

(1)列车为有质量的质点，非匀速运动时需考虑牵引回转系数。

(2)车站抽象为节点，区间线路抽象为弧。

(3)列车可以输出当前速度下不大于最大牵引力的任意大小的牵引力。

文中主要参数定义及变量说明见表1。

表1 主要参数定义及变量说明

2.2 模型及约束条件

目标函数为能耗和时间总体最优，即

minf(X)=(fe(X),ft(X))

(1)

式中：f(X)为总的优化目标；fe(X)为列车消耗的能量；ft(X)为运行所需的时间(X的求解方法将在下一节展开)。这里以s为自变量建立模型。其中，

(2)

受到以下条件约束：

(3)

Fr=μ+κ·v+γ·v2

(4)

Fg≈M·g·θ

(5)

(6)

式中：μ,κ和γ为戴维斯公式中的相关参数；g为重力加速度；θ为线路坡度；A为一个经验常数；R为线路的半径；Ltrain和Lline分别为列车和曲线的长度。

约束式(3)～式(6)为列车运行过程中的受力约束。

(7)

(8)

3 模型求解

考虑到实际线路距离长且线型条件变化频繁，直接划分的计算区段数目太多将消耗较大的计算资源，因此提出一种计算区段简化方法。同时针对基本NSGA-Ⅱ算法收敛速度慢、解分布不够均匀的问题，设计一种以计算区段目标速度为自变量的改进NSGA-Ⅱ算法，并基于变步长仿真获得解所对应的能耗与时间值。

3.1 计算区段简化

由式(4)～式(6)可知，随着列车速度的增加，Fr成为影响运行的主要阻力，因此本文仅对列车速度较大的计算区段进行简化。

(3)加算坡度。针对{p,…,k,…,q}，根据TB/T 1407—1998《列车牵引计算规程》[11]，将曲线附加阻力换算为折算坡度ic，合并相近坡度得化简坡度ig，则加算坡度iall=ic+ig。

(4)根据iall重新划分计算区段。从计算区段p开始，把相同iall的线路划分为一个计算区段，iall值发生改变，则进入下一个计算区段，依此类推，直到计算区段q停止。

3.2 操纵策略与仿真

(3)如果与制动曲线相交，则根据制动曲线操纵。

为了加快仿真速度，采用100、10、1 m变步长仿真。计算距离数组tag，如果区间距离为12 345 m，则tag=[123,4,5]。

Step1初始化：令步长ds=100，累计距离l=0，索引i=0，(s,t,E)=(0, 0, 0)；输入速度v，距离数组tag。

Step2Do:

Fori=0 to 2

Do Whiles

s=s+ds

v1=(v2+ 2 ·a·ds)0.5

t=t+ (v1-v) /a

v=v1

E=E+Ff·ds

EndIf

tag[i+1]=10

s=s-ds

break

EndIf

break

EndIf

EndWhile

ds=ds/10

l=s

EndFor

Step3Ifl

Step4以速度v运行剩下的距离，转Step 5。

Step5输出(v,s,t,E)

其中，Fall为附加阻力。

对于惰行过程

a=-(Fr+Fall)/(1+ρ)·M

(9)

制动过程

a=-(Fb+Fr+Fall)/(1+ρ)·M

(10)

由此即可获得X操纵策略下列车消耗的能量与时间值。

3.3 改进NSGA-Ⅱ算法

将差分进化算法作为交叉算子，同时设计一种三点拥挤距离算子，通过改进NSGA-Ⅱ算法迭代获得能耗-时间Pareto最优的X集合。

(11)

式中：β是控制参数。

判断解的优劣包括两步：①通过非支配排序计算支配等级，等级越低越优；②同等级，则计算拥挤距离为D。

(12)

4 算例研究

以京沪线北京南—廊坊、廊坊—天津南2个区间进行算例研究。线路、车辆等相关数据由项目调研得到。同时，根据现场经验，在计算牵引力时，乘一个牵引力使用系数，见表2。使用Python工具，在64位Windows操作系统上实现高速铁路列车节能驾驶曲线集的生成。

表2 相关参数与值

4.1 线路简化有效性证明

通过化简，北京南—廊坊段计算区段从88个(见图2(a))减少为33个(见图2(b))，仅占原区段数目的37.5%。其中简化后最长计算区段5.10 km，最短0.85 km，平均1.79 km。需要说明的是，因为部分曲线值远大于限速与坡度值，图2(a)中曲线值除以200以便显示。

图2 简化前后线路条件对比

最短时间条件下简化前后的驾驶曲线对比见图3。简化前采用1 m的步长精确仿真，可以看出，2条曲线基本重合。在20 km附近，简化后蓝色实线略滞后于红色实线，是因为相近的线路被合并，而限速取其中较小的值，所以蓝色线保持270 km/h限速运行一段距离，才被允许加速。简化前能耗9.29×109J，简化后为9.28×109J，误差仅占简化前的0.07%。简化前计算时间为23.32 s，简化后时间为1.52 s，计算速度约为原计算方式的15倍。

图3 简化前后最短时间条件下速度-距离曲线

在最短时间条件下运行时，列车能耗最高，但是可以最快地弥补延误。下面选择一种时间较为充足，列车可以采用惰行操纵的情况，进一步对比简化前后仿真的结果。

简化前后T=1 220 S条件下速度-距离曲线见图4。由图4可见，2条曲线整体趋势基本相同，在10～20 km范围偏差较明显，这是因为简化前后的线型条件对应了不同的最优操纵策略。简化前后能耗分别为4.17×109、4.18×109J，相对误差为0.09%；计算时间分别为15.06、1.36 s，计算速度是原计算方式的11倍。

图4 简化前后T=1 220 s条件下速度-距离曲线

通过两种时间条件下仿真结果的对比可知，简化后能耗误差较小，而计算效率大大提升，说明了简化方法的有效性。下面的计算都将在简化线路的基础上进行。

4.2 单区间Pareto曲线

分别用基本NSGA-Ⅱ算法和改进NSGA-Ⅱ算法，计算北京南-廊坊段的Pareto曲线，见图5。图5(a)为基本NSGA-Ⅱ算法所得的不同迭代次数条件下的能耗-时间Pareto曲线，其中第30代到50代产生了较多重复的解，第50代仍未完全收敛且解的分布非常不均匀。图5(b)为改进NSGA-Ⅱ算法的Pareto曲线，从第30代开始，基本趋势已经比较明显，至第50代新生成的解仍在趋势内，因此，我们认为此时已完全收敛。相比传统算法，改进NSGA-Ⅱ的收敛速度更快、解的分布更均匀。

图5 两种方法不同迭代次数的Pareto曲线

为了进一步分析能耗-时间的数值关系，利用matlab中的fit函数拟合所获得的Pareto曲线。Pareto曲线拟合见图6。

图6 Pareto曲线拟合

综合比较，六阶拟合效果最佳，曲线拟合方程为

f1=p1·x6+p2·x5+p3·x4+p4·x3+

p5·x2+p6·x+p7

(13)

式中：p1、p2、p3、p4、p5、p6、p7为系数。

曲线拟合方程系数见表3。

表3 北京南—廊坊段曲线拟合方程系数

同理，获得廊坊—天津南段的Pareto曲线，五阶拟合效果最佳，曲线拟合方程为

f2=p1·x5+p2·x4+p3·x3+

p4·x2+p5·x+p6

(14)

曲线拟合方程系数见表4。

表4 廊坊—天津南段曲线拟合方程系数

4.3 多区间扰动情况下节能驾驶

假设列车从北京南站始发，经停廊坊站，终到天津南站。北京南—廊坊段最小运行时间、计划运行时间、最大运行时间分别为16、21、24 min；廊坊—天津南段分别为15、18、24 min。假设列车因扰动始发晚点5 min。现采用两种策略：①在北京南—廊坊段以最高速度弥补晚点；②保证终到时间不变，利用Pareto曲线寻找可替代节能策略。两种策略条件下的速度-距离曲线见图7。

图7 两种策略条件下的速度-距离曲线

策略1：该策略条件下总能耗为1.78×1010J，总运行时间为34 min。因未获得列车在实际运行条件下的速度-距离曲线，故廊坊—天津南段采用本文提出的对应时间条件下的最优驾驶曲线。

策略2：根据获得的能耗-时间拟合方程，建立数学模型为

minz=f1(x1)+f2(x2)

s.t.x1+x2=34 16≤x1≤24 15≤x2≤24

(15)

通过求解获得，当x1=1 013 s，x2=1 027 s时，线路总能耗为1.66×1010J。

对比两种策略可得：①策略2比策略1节能7.1%；②策略2中列车在2个区间速度分布更加均衡；③两种策略总的运行时间相同，且到达终到站时间均无延误。

5 结论

针对高速铁路线路距离长、线型条件变化频繁，以及列车运行易受扰动偏离基本计划的特点，本文提出了一种高速铁路列车节能驾驶曲线集生成方法，通过简化复杂线路条件，在保证结果准确性的同时提高了计算速度；设计了一种多目标求解算法，基于变步长仿真获得解所对应的目标值；生成一组能耗-时间平衡的Pareto解集，既可以面向常规条件下的列车节能操纵的优化，也可以处理特殊条件下的列车节能调度。通过以上工作得出以下结论：

(1)本文采用的简化频繁变化线路条件的方法，计算区段缩减为原来的37.5%。

(2)以各计算区段目标速度为决策变量，运用变步长仿真法，实例显示，误差仅占简化前逐步精确仿真的0.07%～0.09%，而计算速度是原计算方式的11～15倍。

(3)本文提出的结合差分进化算法和新拥挤距离算子的改进NSGA-Ⅱ算法，收敛速度更快，且解的分布性更好，实例中在第30代收敛趋势明显，第50代已完全收敛。

(4)通过fit函数获得不同区间Pareto前沿的拟合曲线，在扰动发生条件下，较最高速度弥补晚点节能7.1%，为铁路管理者提供了不同运行条件下的节能备选方案，也为节能时刻表的优化奠定了基础。