随机交通网络约束最可靠路径乘子交替方向法

Constrained most reliable path in stochastic traffic network alternating direction method of multipliers

摘要

为分析行程时间相关性和资源约束对交通网络最可靠路径选择的影响，建立考虑行程时间相关性的随机交通网络约束最可靠路径模型，提出了改进乘子交替方向法的求解方法。通过Cholesky分解和二次项放缩，将目标函数转化为可分离结构的矩阵形式，采用块坐标下降法将增广拉格朗日松弛问题分解为一系列子问题，上下界迭代获得近似最优解。针对Sioux Falls网络和Chicago sketch网络进行数值模拟试验并对比分析改进交替方向法和拉格朗日松弛算法的性能。结果表明：有无路段行程时间相关性的约束最可靠路径问题是有区别的；行程时间相关性和资源约束对最可靠路径的选择有显著影响；改进乘子交替方向法的计算效率和收敛性优于拉格朗日松弛算法。

Abstract

In order to analyze the influence of travel time correlation and resource constraints on the selection of the most reliable path in traffic network， the model of the constrained most reliable path with stochastic traffic network considering travel time correlation is established， and the solution method of improved alternating direction method of multipliers is proposed. By Cholesky decomposition and quadratic term deflation， the objective function is transformed into a matrix form of separable structure， and the augmented Lagrangian relaxation problem is decomposed into a series of subproblems by the block coordinate descent method， the upper and lower bounds are iterated to obtain the approximate optimal solution. Numerical simulation experiments are carried out for Sioux Falls network and Chicago sketch network， and the performance of improved alternating direction method of multipliers and Lagrangian relaxation algorithm is compared and analyzed. The results show that the most reliable path problem with or without link travel time correlation is different， the travel time correlation and resource constraints have a significant impact on the selection of the most reliable path， the computational efficiency and convergence of the improved alternating direction method of multipliers are superior to the Lagrangian relaxation algorithm.

Graphical abstract

关键词

交通运输系统工程 / 随机网络 / 约束最可靠路径 / 行程时间相关性 / 乘子交替方向法

Key words

engineering of communication and transportation system / stochastic network / constrained most reliable path / travel time correlation / alternating direction method of multipliers

0 引 言

最优路径问题作为交通网络和车辆导航系统的重要组成部分，在实际交通网络中出行时间是随机的，而且路段之间存在相互作用，行程时间会因交通拥堵、交通事故和天气等而发生关联性变化；同时，路径选择受到行车距离、机动车油耗和电动汽车电量等资源约束的影响。因此，在资源约束和路段行程时间相关条件下，需要对随机交通网络最可靠路径问题进行深入研究，寻找适用于该场景的有效算法。

针对考虑行程时间相关性或资源约束的随机交通网络可靠路径问题已展开研究。Chen等^［1］研究随机时变网络中考虑路段时空相关性的可靠最短路径，以北京市路网为例分析了空间距离、时间距离和道路类型对路段相关性的影响；Rajabi等^［2］通过相关性矩阵量化路段之间的相关性，研究带有时间窗的行程时间相关性的车辆路径问题；Shen等^［3］在信号交通网络中研究考虑路段行程时间相关性，交叉口和交通信号延迟相关性的最可靠路径问题；Zockaie等^［4］以准时到达概率和最小行程时间为可靠性标准，在动态随机网络中考虑行程时间的相关性；Zhang等^［5］基于高速公路网络和城市道路网络的真实数据，研究了具有时间相关性的随机最短路问题；Pugliese等^［6］研究了约束最可靠路径的变化形式，其中资源消耗是不确定参数，展示可变性对最可靠路径选择的影响；Nie等^［7］针对燃料消耗和温室气体排放构建了约束最短路模型，考虑了不同约束对最优路径的影响；潘义勇等^［8］以期望-方差为可靠性标准，研究资源约束条件下的路径选择行为；Zeng等^［9］将考虑准时到达概率和燃油消耗作为约束条件，寻找可靠且具有最优油耗的路径；Thomas等^［10］采用双向A*算法在多种大型交通网络中求解带资源约束的最短路径问题。上述研究集中在行程时间相关性或约束最可靠路径问题上，但同时考虑这两个方面更符合现实交通网络。因此，有必要进一步探索考虑行程时间相关性的随机交通网络约束最可靠路径问题。

目前求解最优路径问题主要采用拉格朗日松弛算法、剪枝算法、A*算法等。Prakash^［11］提出了基于边界的双层剪枝准则，结合标签校正算法确定随机时变网络上的最可靠路径；Yang等^［12］以准时到达概率和百分位旅行时间为可靠性标准，拓展拉格朗日分解方法，重构了可靠路径模型。拉格朗日松弛方法将复杂的约束松弛到目标函数中，使复杂问题转化为更易处理的子问题，适用于具有复杂约束的优化问题，但在处理大规模优化问题时，会出现计算量过大的问题，导致收敛速度降低；剪枝算法可以排除非优路径，减少需要搜寻的路径，适用于动态规划中求解最优路径，其效果依赖于剪枝策略的选择，当最优解的部分路径被删除时导致解的质量下降或错过最优解，其解与最优解之间的差距也是未知的；乘子交替方向法（Alternating direction method of multiploers，ADMM）结合增广拉格朗日松弛法和块坐标下降法，将分解的子问题并行求解，在处理大规模优化问题时具有较高的计算效率，同时保证了良好的收敛性。目前ADMM方法在交通领域应用较少且主要集中在VRP问题和车辆充电桩选址问题，Song等^［13］基于增广拉格朗日松弛方法，通过两次拉格朗日松弛构造均值-标准差的车辆路径问题的下界，得到高质量的可行解；Yao等^［14］利用ADMM方法分解带时间窗的VRP问题，交替更新乘子和求解每个子问题，提高原始解和对偶解的质量；Song等^［15］将ADMM方法引入可靠路径问题中，基于ADMM方法分解目标函数，提高了解的收敛性和解的质量。而ADMM方法还未用于求解带有行程时间相关性的约束最可靠路径，基于ADMM构建新的分解方法求解该问题，可对时间相关网络下约束最可靠路径算法进行拓展。

综上所述，本文同时考虑行程时间相关性和资源约束条件，建立考虑行程时间相关性的随机交通网络约束最可靠路径模型，提出了基于乘子交替方向法的分解方法，构建相应的求解算法。然而，由于二次罚项和协方差相关系数矩阵导致了目标函数的非线性，难以直接分解，这也是算法构建的重点。第一，定义可靠性标准以及协方差相关系数作为路径目标函数，建立考虑行程时间相关性的约束最可靠路径模型；第二，采用Cholesky分解和二次项放缩将增广拉格朗日松弛问题分解为子问题；第三，求解子问题，得到最优路径解的上界和下界，并通过路径解迭代和乘子更新不断逼近最优解；第四，编写Python代码算法程序，针对交通网络Sioux Falls network和Chicago sketch network进行数值试验，验证算法可行性，对比分析基于ADMM方法和基于LR算法的性能；第五，总结和未来可能的研究方向。

1 问题描述与建模

给定交通网络拓扑图G=\left(N,A,X\right)，N=\left(\left|N\right|=n\right)表示交通节点的集合，A=\left(\left|A\right|=m\right)表示路段的集合，X表示路段权值的集合，每条路段ij对应一个行程时间t_{ij}和资源消耗w_{ij}，其中行程时间t_{ij}是随机的，其期望值c_{ij}、标准差{\sigma }_{ij}和表示路段ij和kl的行程时间协方差\mathrm{C}\mathrm{o}\mathrm{v}\left(ij,kl\right)可通过路网数据获得，以期望-均方差作为可靠性标准，考虑路段行程时间相关性的随机网络约束最可靠路径数学模型如下：

\mathrm{M}\mathrm{i}\mathrm{n}\sum _{ij\in A}{c}_{ij}{x}_{ij}+\lambda \sqrt[]{\sum _{ij\in A}{\sigma }_{ij}^{\mathrm{2}}{x}_{ij}+\sum _{ij,kl\in A}\mathrm{C}\mathrm{o}\mathrm{v}\left(ij,kl\right)x_{ij}{x}_{kl}}（1）

\sum _{(i,j\in A)}{w}_{ij}{x}_{ij}\le W（2）

\mathrm{s}.\mathrm{t}.\sum _{j:ij\in A}{x}_{ij}-\sum _{j:ji\in A}{x}_{ij}=\left\{\begin{array}{c}\mathrm{1},i=o\\ \mathrm{0},i\in N-\left\{o,d\right\}\\ -\mathrm{1},i=d\end{array}\right.（3）

x_{ij}\in \left\{\mathrm{0,1}\right\}（4）

其中，目标函数（1）考虑行程时间的可变性和行程时间相关性，使所有车辆的总路线成本最低；约束（2）确保车辆资源消耗限制，w_{ij}为路段ij的资源消耗，W为路径资源总消耗的上限值，如汽车燃油消耗总量和电池剩余电量等；约束（3）保证网络的流量平衡；约束（4）保证决策变量是二元的，x_{ij}=\mathrm{1}表示路段ij为最可靠路径的组成部分，\lambda为可靠性系数。

依据x_{ij}的二元变量特性，可将目标函数转化为矩阵形式，得到约束最可靠路径的矩阵形式如下：

\mathrm{M}\mathrm{i}\mathrm{n}\hspace{1em}{\mathit{c}}^{\mathrm{T}}\mathit{x}+\lambda \sqrt[]{{\mathit{x}}^{\mathrm{T}}\mathit{\Sigma }\mathit{x}}（5）

\mathrm{s}.\mathrm{t}.\left\{\begin{array}{l}{\mathit{w}}^{\mathrm{T}}\mathit{x}\le W\\ \mathrm{式}\left(\mathrm{3}\right),\mathrm{式}\left(\mathrm{4}\right)\end{array}\right.（6）

式中：矩阵\mathit{c}={\left\{\begin{array}{cccc}{c}_{ij}^{\mathrm{1}}& c_{ij}^{\mathrm{2}}& \cdots & c_{ij}^n\end{array}\right\}}^{\mathrm{T}}为路段的期望值；矩阵\mathit{w}={\left\{\begin{array}{cccc}{w}_{ij}^{\mathrm{1}}& w_{ij}^{\mathrm{2}}& \cdots & w_{ij}^n\end{array}\right\}}^{\mathrm{T}}的元素为路段的资源消耗； \mathit{\Sigma }=\left\{\begin{array}{cccc}{\sigma }_{\mathrm{11}}^{\mathrm{2}}& {\rho }_{\mathrm{12}}& \cdots & {\rho }_{\mathrm{1}n}\\ {\rho }_{\mathrm{21}}& {\sigma }_{\mathrm{22}}^{\mathrm{2}}& \cdots & {\rho }_{\mathrm{2}n}\\ ⋮& ⋮& \ddots & ⋮\\ {\rho }_{n\mathrm{1}}& {\rho }_{n\mathrm{2}}& \cdots & {\sigma }_{nn}^{\mathrm{2}}\end{array}\right\}为方差-协方差矩阵，主对角线的元素为路段行程时间的方差，其他元素为路段ij和kl的协方差。

\mathit{\Sigma }具有半正定矩阵的性质，可以通过Cholesky分解为下三角矩阵\mathit{P}和上三角矩阵，分解形式如下：

\mathit{\Sigma }=\mathit{P}{\mathit{P}}^{\mathrm{T}}=\left\{\begin{array}{cccc}{p}_{\mathrm{11}}& \mathrm{0}& \cdots & \mathrm{0}\\ p_{\mathrm{21}}& p_{\mathrm{22}}& \cdots & \mathrm{0}\\ ⋮& ⋮& \ddots & ⋮\\ p_{n\mathrm{1}}& p_{n\mathrm{2}}& \cdots & p_{nn}\end{array}\right\}\left\{\begin{array}{cccc}{p}_{\mathrm{11}}& p_{\mathrm{21}}& \cdots & p_{n\mathrm{1}}\\ \mathrm{0}& p_{\mathrm{22}}& \cdots & p_{n\mathrm{2}}\\ ⋮& ⋮& \ddots & ⋮\\ \mathrm{0}& \mathrm{0}& \cdots & p_{nn}\end{array}\right\}（7）

目标函数中根号项与协方差的存在导致难以分解，引入新的变量\mathit{y}和等式约束\mathit{y}={\mathit{P}}^{\mathrm{T}}\mathit{x}，对目标函数进行重构，得到约束最可靠路径问题的等价形式为：

\mathrm{M}\mathrm{i}\mathrm{n}{\mathit{c}}^{\mathrm{T}}\mathit{x}+\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}（8）

\mathrm{s}.\mathrm{t}.\left\{\begin{array}{l}{\mathit{w}}^{\mathrm{T}}\mathit{x}\le W\\ \mathit{y}={\mathit{P}}^{\mathrm{T}}\mathit{x}\\ \mathrm{式}\left(\mathrm{3}\right),\mathrm{式}\left(\mathrm{4}\right)\end{array}\right.（9）

2 乘子交替方向法分解

2.1　拉格朗日松弛

将式（9）的约束条件进行拉格朗日松弛，引入拉格朗日乘子\mathit{\mu }和v，可以得到标准的拉格朗日松弛问题：

\begin{array}{c}{L}_{\rho }\left(\mathit{x},\mathit{y},\mathit{\mu },v\right)=\mathrm{M}\mathrm{i}\mathrm{n}{\mathit{c}}^{\mathrm{T}}\mathit{x}+\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}+\\ v\left({\mathit{w}}^{\mathrm{T}}\mathit{x}-W\right)+{\mathit{\mu }}^{\mathrm{T}}\left({\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathit{y}\right)\end{array}（10）

对其变量重组与分解，得到关于\mathit{x} 和\mathit{y}的两个子问题：

\begin{array}{c}L\left(\mathit{x},\mathit{y},\mathit{\mu },v\right)=\mathrm{M}\mathrm{i}\mathrm{n}{\left(\mathit{c}+\mathit{P}\mathit{\mu }+v\mathit{w}\right)}^{\mathrm{T}}\mathit{x}+\\ \lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}-{\mathit{\mu }}^{\mathrm{T}}\mathit{y}-vW\end{array}（11）

L\left(\mathit{x},\mathit{\mu },v\right)=\underset{x}{\mathrm{M}\mathrm{i}\mathrm{n}}{\left(\mathit{c}+\mathit{P}\mathit{\mu }+v\mathit{w}\right)}^{\mathrm{T}}\mathit{x}（12）

L\left(\mathit{y},\mathit{\mu },v\right)=\underset{y}{\mathrm{M}\mathrm{i}\mathrm{n}}\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}-{\mathit{\mu }}^{\mathrm{T}}\mathit{y}-vW（13）

为了得到拉格朗日松弛问题的最大下界，推导对偶问题如下：

\underset{\mathrm{\mu },\mathrm{\nu }\in R^{\mathrm{m}}}{\mathrm{m}\mathrm{a}\mathrm{x}}L\left(\mathit{x},\mathit{y},\mathit{\mu },v\right)\le l\left(\mathit{\mu },v\right)（14）

通过求解其分解的两个子问题，得到原问题最优值的下界。

迭代过程根据式（15）和式（16）对拉格朗日乘子进行更新：

{\mathit{\mu }}^{k+\mathrm{1}}={\mathit{\mu }}^k+{\theta }_{\mathrm{1}}\left({\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathit{y}\right)（15）

v^{k+\mathrm{1}}=v^k+{\theta }_{\mathrm{2}}\left({\mathit{w}}^{\mathrm{T}}\mathit{x}-W\right)（16）

每次迭代的步长为：

{\theta }_{\mathrm{1}}=\frac{\eta \left(\mathrm{U}\mathrm{B}-\mathrm{L}\mathrm{B}\right)}{{‖{\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathit{y}‖}^{\mathrm{2}}}（17）

{\theta }_{\mathrm{2}}=\frac{\eta \left(\mathrm{U}\mathrm{B}-\mathrm{L}\mathrm{B}\right)}{{‖{\mathit{w}}^{\mathrm{T}}\mathit{x}-W‖}^{\mathrm{2}}}（18）

式中：\mathrm{0}\le \eta \le \mathrm{2}。

2.2　增广拉格朗日松弛

为了得到增广拉格朗日松弛问题，引入二次罚项，\rho \ge \mathrm{0}是二次惩罚项的惩罚参数，构造的增广拉格朗日松弛问题如下：

\begin{array}{l}{L}_{\rho }\left(\mathit{x},\mathit{y},\mathit{\mu },v\right)=\mathrm{M}\mathrm{i}\mathrm{n}{\mathit{c}}^{\mathrm{T}}\mathit{x}+\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}+v\left({\mathit{w}}^{\mathrm{T}}\mathit{x}-W\right)+\\ {\mathit{\mu }}^{\mathrm{T}}\left({\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathit{y}\right)+\frac{\rho }{\mathrm{2}}{‖{\mathit{w}}^{\mathrm{T}}\mathit{x}-W‖}_{\mathrm{2}}^{\mathrm{2}}+\frac{\rho }{\mathrm{2}}{‖{\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathit{y}‖}_{\mathrm{2}}^{\mathrm{2}}\end{array}（19）

ADMM方法是增广拉格朗日松弛法和块坐标下降法的结合，可用于解决具有两个目标函数的凸优化问题，在基于ADMM分解方法的框架中，{\mathit{x}}^{k+\mathrm{1}}为第k+\mathrm{1}次迭代的路径解，{\mathit{y}}^{k+\mathrm{1}}为第k+\mathrm{1}次迭代的辅助变量，{\mathit{\mu }}^{k+\mathrm{1}}和v^{k+\mathrm{1}}为乘子更新，求解最可靠路径{\mathit{x}}^{k+\mathrm{1}}时变量{\mathit{y}}^k作为固定项，而在求解变量{\mathit{y}}^{k+\mathrm{1}}时{\mathit{x}}^{k+\mathrm{1}}作为固定项，两者交替更新。每次迭代时需要进行两部分优化，子问题解迭代和拉格朗日乘子更新。

{\mathit{x}}^{k+\mathrm{1}}=\underset{x}{\mathrm{a}\mathrm{r}\mathrm{g}\mathrm{m}\mathrm{i}\mathrm{n}}{L}_{\rho }\left(\mathit{x},{\mathit{y}}^k,{\mathit{\mu }}^k,v^k\right)（20）

{\mathit{y}}^{k+\mathrm{1}}=\underset{y}{\mathrm{a}\mathrm{r}\mathrm{g}\mathrm{m}\mathrm{i}\mathrm{n}}{L}_{\rho }\left({\mathit{x}}^{k+\mathrm{1}},\mathit{y},{\mathit{\mu }}^k,v^k\right)（21）

{\mathit{\mu }}^{k+\mathrm{1}}={\mathit{\mu }}^k+\rho \left({\mathit{P}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}-{\mathit{y}}^{k+\mathrm{1}}\right)（22）

v^{k+\mathrm{1}}=v^k+\rho \left({\mathit{w}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}-W\right)（23）

增广拉格朗日具有较好的解可行性和收敛性，但二次罚项存在，无法直接求解最可靠路径，展开二次项如下：

{‖{\mathit{w}}^{\mathrm{T}}\mathit{x}-W‖}_{\mathrm{2}}^{\mathrm{2}}={\left({\mathit{w}}^{\mathrm{T}}\mathit{x}-W\right)}^{\mathrm{2}}={\left({\mathit{w}}^{\mathrm{T}}\mathit{x}\right)}^{\mathrm{2}}-\mathrm{2}W{\mathit{w}}^{\mathrm{T}}\mathit{x}+W^{\mathrm{2}}（24）

{‖{\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathit{y}‖}_{\mathrm{2}}^{\mathrm{2}}={\mathit{x}}^{\mathrm{T}}\mathit{P}{\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathrm{2}{\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}\mathit{x}+{\mathit{y}}^{\mathrm{T}}\mathit{y}（25）

二元变量平方等于它本身，即{\left(x_{ij}\right)}^{\mathrm{2}}=x_{ij}，依据该性质对式（24）和式（25）进行放缩：

\begin{array}{c}{\left({\mathit{w}}^{\mathrm{T}}\mathit{x}\right)}^{\mathrm{2}}={\left(\sum _{ij\in A}{w}_{ij}{x}_{ij}\right)}^{\mathrm{2}}\le \sum _{ij\in A}m{w}_{ij}^{\mathrm{2}}{x}_{ij}^{\mathrm{2}}=\\ \sum _{ij\in A}m{w}_{ij}^{\mathrm{2}}{x}_{ij}=m{\mathit{w}}_{\mathrm{1}}^{\mathrm{T}}\mathit{x}\end{array}（26）

\begin{array}{c}{\mathit{x}}^{\mathrm{T}}\mathit{P}{\mathit{P}}^{\mathrm{T}}\mathit{x}={\sum _{ij\in A}\left(\sum _{kl\in A}{p}_{ij,kl}{x}_{ij}\right)}^{\mathrm{2}}\le \\ m^{\mathrm{2}}\sum _{ij\in A}\sum _{kl\in A}{p}_{ij.kl}^{\mathrm{2}}{x}_{ij}=m^{\mathrm{2}}{\mathit{P}}_{\mathrm{1}}^{\mathrm{T}}\mathit{x}\end{array}（27）

式中：m为网络路段数，由于x_{ij}的二元变量特性，部分w_{ij} 和p_{ij,kl}不会同时取到；矩阵{\mathit{w}}_{\mathrm{1}}的元素为路段资源约束的平方w_{ij}^{\mathrm{2}}；矩阵{\mathit{P}}_{\mathrm{1}}的元素为Cholesky分解矩阵的分量最大值的平方。

从而得到放缩之后的约束最可靠路径增广拉格朗日松弛函数的表达式：

\begin{array}{c}{L}_{\rho }\left(\mathit{x},\mathit{y},\mathit{\mu },v\right)=\mathrm{M}\mathrm{i}\mathrm{n}{\mathit{c}}^{\mathrm{T}}\mathit{x}+\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}+{\mathit{\mu }}^{\mathrm{T}}\left({\mathit{P}}^{\mathrm{T}}\mathit{x}-\mathit{y}\right)+\\ v\left({\mathit{w}}^{\mathrm{T}}\mathit{x}-W\right)+\frac{\rho }{\mathrm{2}}\left(m{\mathit{w}}_{\mathrm{1}}^{\mathrm{T}}\mathit{x}-\mathrm{2}W{\mathit{w}}^{\mathrm{T}}\mathit{x}+W^{\mathrm{2}}\right)+\\ \frac{\rho }{\mathrm{2}}\left(m^{\mathrm{2}}{\mathit{P}}_{\mathrm{1}}^{\mathrm{T}}\mathit{x}-\mathrm{2}{\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}\mathit{x}+{\mathit{y}}^{\mathrm{T}}\mathit{y}\right)\end{array}（28）

对增广拉格朗日松弛函数L_{\rho }\left(\mathit{x},\mathit{y},\mathit{\mu },v\right)进行变量重组与分解，得到关于\mathit{x}和\mathit{y}两个子问题：

\begin{array}{c}\underset{x}{\mathrm{a}\mathrm{r}\mathrm{g}\mathrm{m}\mathrm{i}\mathrm{n}}{L}_{\rho }\left(\mathit{x},{\mathit{y}}^k,{\mathit{\mu }}^k,v^k\right)=(\mathit{c}+\mathit{P}{\mathit{\mu }}^k+v^k\mathit{w}+\\ \frac{\rho }{\mathrm{2}}m{\mathit{w}}_{\mathrm{1}}-\rho W\mathit{w}+\frac{\rho }{\mathrm{2}}{m}^{\mathrm{2}}{\mathit{P}}_{\mathrm{1}}-\rho \mathit{P}{\mathit{y}}^k{)}^{\mathrm{T}}\mathit{x}\end{array}（29）

\begin{array}{c}\underset{y}{\mathrm{a}\mathrm{r}\mathrm{g}\mathrm{m}\mathrm{i}\mathrm{n}}{L}_{\rho }\left({\mathit{x}}^{k+\mathrm{1}},\mathit{y},{\mathit{\mu }}^k,v^k\right)=\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}-\\ {\left({\mathit{\mu }}^k\right)}^{\mathrm{T}}\mathit{y}+{\mathit{y}}^T\mathit{y}-\rho {\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}\end{array}（30）

关于\mathit{x}的子问题为最短路径问题，其中含有\mathit{y}项的\mathit{P}{\mathit{y}}^k作为固定项，可以采用标签更正算法求解。关于\mathit{y}的子问题，参数\lambda、\mathit{\mu }和\rho初始给定，{\mathit{P}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}作为固定项，令{‖\mathit{y}‖}_{\mathrm{2}}=t\ge \mathrm{0}，\mathrm{m}\mathrm{a}\mathrm{x}\left\{\left.{\mathit{\mu }}^{\mathrm{T}}\mathit{y}+\rho {\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}\mathit{x}:{‖\mathit{y}‖}_{\mathrm{2}}=t\right\}\right.={‖\mathit{\mu }+\rho {\mathit{P}}^{\mathrm{T}}\mathit{x}‖}_{\mathrm{2}}t，最优值\mathit{y}=\frac{t\left(\mathit{\mu }+\rho {\mathit{P}}^{\mathrm{T}}\mathit{x}\right)}{‖{\mathit{\mu }+\rho {\mathit{P}}^{\mathrm{T}}\mathit{x}‖}_{\mathrm{2}}}。则：

\begin{array}{c}\mathrm{m}\mathrm{i}\mathrm{n}\left\{\left.\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}-{\mathit{\mu }}^{\mathrm{T}}\mathit{y}+{\mathit{y}}^{\mathrm{T}}\mathit{y}-\rho {\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}\mathit{x}\right\}\right.=\\ \underset{t\ge \mathrm{0}}{\mathrm{m}\mathrm{i}\mathrm{n}}\underset{‖\mathit{y}‖=t}{\mathrm{m}\mathrm{i}\mathrm{n}}\left\{\left.\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}-{\mathit{\mu }}^{\mathrm{T}}\mathit{y}+{\mathit{y}}^{\mathrm{T}}\mathit{y}-\rho {\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}\mathit{x}\right\}=\right.\\ \underset{t\ge \mathrm{0}}{\mathrm{m}\mathrm{i}\mathrm{n}}\left\{\left.\lambda t+t^{\mathrm{2}}-\underset{‖\mathit{y}‖=t}{\mathrm{m}\mathrm{a}\mathrm{x}}\left({\mathit{\mu }}^{\mathrm{T}}\mathit{y}+\rho {\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}\mathit{x}\right)\right\}\right.=\\ \underset{t\ge \mathrm{0}}{\mathrm{m}\mathrm{i}\mathrm{n}}\left\{\left.\lambda t+t^{\mathrm{2}}-{‖\mathit{\mu }+\rho {\mathit{P}}^{\mathrm{T}}\mathit{x}‖}_{\mathrm{2}}t\right\}\right.\end{array}（31）

可得：

\begin{array}{l}\mathrm{m}\mathrm{i}\mathrm{n}\left\{\left.\lambda \sqrt[]{{\mathit{y}}^{\mathrm{T}}\mathit{y}}-{\left({\mathit{\mu }}^k\right)}^{\mathrm{T}}\mathit{y}+{\mathit{y}}^{\mathrm{T}}\mathit{y}-\rho {\mathit{y}}^{\mathrm{T}}{\mathit{P}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}\right\}\right.=\\ \left\{\begin{array}{l}\mathrm{0},{‖{\mathit{\mu }}^k+\rho {\mathit{P}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}‖}_{\mathrm{2}}\le \lambda \\ -\frac{{\left({‖{\mathit{\mu }}^k+\rho {\mathit{P}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}‖}_{\mathrm{2}}-\lambda \right)}^{\mathrm{2}}}{\mathrm{4}}, {‖{\mathit{\mu }}^k+\rho {\mathit{P}}^{\mathrm{T}}{\mathit{x}}^{k+\mathrm{1}}‖}_{\mathrm{2}}\ge \lambda \end{array}\right.\end{array}（32）

3 求解算法

本文采用改进乘子交替方向法求解考虑行程时间相关性的随机交通网络约束最可靠路径问题。首先进行参数初始化，包括上界和下界、拉格朗日乘子和惩罚参数；其次求解子问题式（29）得到路径解{\mathit{x}}^{k+\mathrm{1}}，再基于路径解{\mathit{x}}^{k+\mathrm{1}}求解子问题式（30）得到y^{k+\mathrm{1}}；之后，将解{\mathit{x}}^{k+\mathrm{1}}代入式（5）得到上界，以{\mathit{\mu }}^k和v^k为基础求解标准拉格朗日松弛问题式（10）获得下界；最后，更新拉格朗日乘子和惩罚参数。

基于改进乘子交替方向法分解的可靠路径求解算法流程如下：

步骤 1 初始化

设置迭代次数k=\mathrm{1}，最大迭代次数k_{\mathrm{m}\mathrm{a}\mathrm{x}}

设置拉格朗日乘子\mathit{\mu }和v以及惩罚参数\rho

设置上界\mathrm{U}\mathrm{B}=\mathrm{100}和下界\mathrm{L}\mathrm{B}=-\mathrm{50}，

步骤 2 求解可行解，计算上下界

求解L_{\rho }\left(\mathit{x},{\mathit{y}}^k,{\mathit{\mu }}^k,v^k\right)得到解{\mathit{x}}^{k+\mathrm{1}}

基于路径解{\mathit{x}}^{k+\mathrm{1}}，求解L_{\rho }\left({\mathit{x}}^{k+\mathrm{1}},\mathit{y},{\mathit{\mu }}^k,v^k\right)得到{\mathit{y}}^{k+\mathrm{1}}

求解标准的拉格朗日松弛问题式（10）得到下界\mathrm{L}\mathrm{B}

将{\mathit{x}}^{k+\mathrm{1}}代入式（5）中，计算上界\mathrm{U}\mathrm{B}

步骤 3 拉格朗日乘子更新

更新乘子{\mathit{\mu }}^{k+\mathrm{1}}和v^{k+\mathrm{1}}

路径解\mathit{x}、变量\mathit{y}、下界\mathrm{L}\mathrm{B}和上界\mathrm{U}\mathrm{B}迭代

步骤 4 算法终止\frac{\mathrm{U}\mathrm{B}-\mathrm{L}\mathrm{B}}{\mathrm{U}\mathrm{B}}\le \widehat{\epsilon }

if\frac{\mathrm{U}\mathrm{B}-\mathrm{L}\mathrm{B}}{\mathrm{U}\mathrm{B}}\le \widehat{\epsilon } then

终止算法

输出\mathit{x}和误差\epsilon

else

k=k+\mathrm{1}，返回步骤 2

4 数值实验

为验证上述模型和算法的有效性和性能，采用小型测试网络Sioux Falls network和大网络Chicago sketch network进行数值实验。Sioux Falls network包含24个节点和76条边，网络规模较小且网络数据易获取，是标准的交通测试网络，适用于对基于ADMM的方法进行初步验证；Chicago sketch network基于真实的城市交通网络，包含933个节点和2 950条边，代表大规模交通网络，适用于对比分析基于ADMM方法和LR算法在大规模网络上的性能。为模拟交通网络的随机特性，构造带有资源约束和路段行程时间相关性的随机交通网络，对路段的相关权值进行数值模拟，路段行程时间的期望c_{ij}和方差{\sigma }_{ij}^{\mathrm{2}}在一定区间内随机取得，路段的资源约束w_{ij}从（0，20］随机取得，如表1所示，两种方法的最大迭代次数为50次，初始值上界UB=100、下界LB=-50，可靠性系数\lambda =\mathrm{2}，相对间隙设为5.00%。

4.1　小网络

首先，对模型中的可靠性系数和惩罚参数敏感度进行分析，如表2所示。\lambda =\mathrm{2}和\rho =\mathrm{0.1}时上界UB=37.42，迭代次数为38次；\lambda =\mathrm{0.5}和\rho =\mathrm{0.1}时上界UB=26.03，迭代次数为11次；\lambda =\mathrm{2}和\rho =\mathrm{1}时上界UB=37.42，迭代次数为5次；可靠性系数\lambda为路径规划中的风险偏好，当\lambda较大时，需寻找更可靠的路径，导致更多的上界和迭代次数；惩罚参数\rho显著影响迭代次数和收敛速度，当\rho过大时，会出现迭代过度和收敛过快的情况，错过其他路径解，导致解的质量下降。

在小型网络Sioux Falls network（见图1）中验证基于ADMM的分解方法的可行性。图2是起终点为1→19的求解过程中上下界迭代变化图，随着逐步迭代，上界UB在第3次减小后保持不变，代表在迭代的过程中最优路径解变化，下界LB逐渐增大，上下界相对间隙减小不断逼近最优解，最终相对间隙小于5.00%。这符合考虑路段行程时间相关性的约束最可靠路径模型及其算法的理论思路，同时验证了算法的收敛性，证明本文推导的改进乘子交替方向法的分解方法是正确的，该算法可用于求解具有路段行程时间相关性的随机交通网络约束最可靠路径问题。

对路段行程时间相关性条件下不同资源约束上限的可靠路径问题和相同资源约束上限情况下有无相关性的可靠路径问题进行求解分析。表3选取3个相同起终点在不同资源约束上限W的情况下，最可靠路径的选择、最优值的上下界以及相对间隙，其中W=+\mathrm{\infty }表示不考虑资源约束上限时路段行程时间相关的最可靠路径模型求解。由表3可得，起终点为1→23在不同的资源约束上限W=+\mathrm{\infty }和W=\mathrm{29}的情况下，求解的最可靠路径是不同的；起终点为3→17在不同的资源约束上限W=+\mathrm{\infty }和W=\mathrm{28}的情况下，求解的最可靠路径是不同的；其他路径在不同资源约束上限的情况下也有不同的最可靠路径。

表4为相同起终点和相同资源约束上限时是否具有路段行程时间相关性的最可靠最短路径选择对比，起终点为3→17在资源约束上限W=+\mathrm{\infty }的条件下，具有路段行程时间相关性和无路段行程时间相关性的最可靠路径是不同的；起终点为4→18在资源约束上限W=+\mathrm{\infty }的条件下，具有路段行程时间相关性和无路段行程时间相关性的最可靠路径是不同的；起终点为4→18在资源约束上限W=\mathrm{25}的条件下，具有路段行程时间相关性和无路段行程时间相关性的最可靠路径是不同的，考虑路段行程时间相关性对随机交通网络的最可靠路径求解有影响。

4.2　大网络

芝加哥网络如图3所示，采用基于ADMM的方法在其网络中求解约束最可靠路径，并与基于拉格朗日松弛（LR）算法的性能对比分析。选取两个不同路径进行求解，其上界、上界和相对间隙的变化如图4和图5所示。图4中迭代到第16次时，基于ADMM的方法相对间隙小于5.0%，而基于LR算法的相对间隙大于5.0%，继续迭代到第25次时，最终相对间隙小于5.0%，两者相比ADMM的方法能更快地确定约束最可靠路径的最优解且收敛速度更快。图5中两种方法迭代30次，基于ADMM的方法上界小于基于LR算法的上界，基于ADMM的方法相对间隙为1.63%，而基于LR算法的相对间隙为4.22%，前者的相对间隙收敛速度较快，收敛结果优于基于LR算法。从生成的上界、下界变化速度以及相对间隙的变化来看，基于ADMM的方法优于LR算法，收敛速度更快。

选取不同的路径进行CPU运行时间、迭代次数和上下界对比，如表5所示，基于ADMM方法的CPU运行时间总体上少于基于LR算法，如路径1和路径3的测试中相同的迭代次数，基于ADMM方法的CPU运行时间较少，分别减少了14.6%和9.6%，且具有更小的相对间隙；路径5的测试中，基于LR算法的CPU运行时间和迭代次数分别为67 s和33次，基于ADMM方法的CPU运行时间和迭代次数分别为59 s和29次，CPU运行时间减少了11.9%，需要的迭代次数更少。基于ADMM方法和基于LR算法的上界通过式（5）确定，基于ADMM方法的上界总体上小于基于算法LR的上界，其最终解更趋向于最优解。从CPU运行时间、迭代次数和上界值来看，多数场景下基于ADMM的方法性能优于LR算法，能够快速获取高质量可行解并保持更好的收敛性。

由此可得，基于ADMM的分解方法适用于求解具有路段行程时间相关性的随机交通网络约束最可靠路径问题；考虑路段行程时间相关性和不考虑路段行程时间相关性的约束最可靠路径问题是有根本区别的，路段行程时间相关性和资源约束对路径的选择有显著影响。在大型网络中对比了基于ADMM的方法和基于LR算法的上下界和相对间隙的变化，同时对比了CPU运行时间和迭代次数，基于ADMM的分解方法计算效率和收敛性优于LR算法，能够快速获得可行解。

5 结 论

（1）本文建立的最可靠路径模型能同时反映交通网络中路段之间相关特性和资源约束的影响，更符合实际交通网络路径选择情形。数值实验表明：有无路段行程时间相关性的约束最可靠路径问题是有区别的；行程时间相关性和资源约束对最可靠路径的选择有显著影响。

（2）提出的基于ADMM的分解方法能较快获取高质量可行解并保持更好的收敛性。从生成的上下界、相对间隙的变化以及CPU运行时间和迭代次数来看，基于ADMM的分解方法其计算效率和收敛性优于拉格朗日松弛算法。

参考文献

[1] Chen P, Tong R, Yu B, et al. Reliable shortest path finding in stochastic time-dependent road network with spatial-temporal link correlations: a case study from Beijing[J]. Expert Systems with Applications, 2020, 147: 113192.

[2] Rajabi B M, Shariat M A, Babaei M, et al. Reliable vehicle routing problem in stochastic networks with correlated travel times[J]. Operational Research, 2021, 21: 299-330.

[3] Shen L, Shao H, Wu T, et al. Finding the reliable shortest path with correlated link travel times in signalized traffic networks under uncertainty[J]. Transportation Research Part E: Logistics and Transportation Review, 2020, 144: 02159.

[4] Zockaie A, Mahmassani H S, Nie Y. Path finding in stochastic time varying networks with spatial and temporal correlations for heterogeneous travelers[J]. Transportation Research Record, 2016(1): 105-113.

[5] Zhang D, Wallace S W, Guo Z, et al. On scenario construction for stochastic shortest path problems in real road networks[J]. Transportation Research Part E: Logistics and Transportation Review, 2021, 152: 102410.

[6] Pugliese L D P, Guerriero F, Poss M. The resource constrained shortest path problem with uncertain data: a robust formulation and optimal solution approach[J]. Computers & Operations Research, 2019, 107: 140-155.

[7] Nie Y M, Li Q. An eco-routing model considering microscopic vehicle operating conditions[J]. Transportation Research Part B: Methodological, 2013, 55: 154-170.

[8] 潘义勇, 马健霄. 基于可靠性的随机交通网络约束最优路径问题[J]. 东南大学学报:自然科学版, 2017, 47(6): 1263-1268.

[9] Pan Yi-yong, Ma Jian-xiao. Constrained shortest path problem in stochastic traffic network based on reliability[J]. Journal of Southeast University (Natural Science Edition), 2017, 47(6): 1263-1268.

[10] Zeng W, Miwa T, Morikawa T. Eco-routing problem considering fuel consumption and probabilistic travel time budget[J]. Transportation Research Part D: Transport and Environment, 2020, 78: 102219.

[11] Thomas B W, Calogiuri T, Hewitt M. An exact bidirectional A* approach for solving resource-constrained shortest path problems[J]. Networks, 2019, 73(2): 187-205.

[12] Prakash A A. Algorithms for most reliable routes on stochastic and time-dependent networks[J]. Transportation Research Part B: Methodological, 2020, 138: 202-220.

[13] Yang L, Zhou X. Optimizing on-time arrival probability and percentile travel time for elementary path finding in time-dependent transportation networks: linear mixed integer programming reformulations[J]. Transportation Research Part B: Methodological, 2017, 96: 68-91.

[14] Song M, Cheng L. An augmented lagrangian relaxation method for the mean-standard deviation based vehicle routing problem[J]. Knowledge-Based Systems, 2022, 247: 108736.

[15] Yao Y, Zhu X, Dong H, et al. ADMM-based problem decomposition scheme for vehicle routing problem with time windows[J]. Transportation Research Part B: Methodological, 2019, 129: 156-174.

[16] Song M, Cheng L. Solving the reliability-oriented generalized assignment problem by lagrangian relaxation and alternating direction method of multipliers[J]. Expert Systems with Applications, 2022, 205: 117644.

基金资助

国家自然科学基金项目(51508280)