An Enhanced Ant Colony Optimization Combined With Clustering Decomposition for Solving Complex Green Vehicle Routing Problem
-
摘要: 针对带时间窗的低能耗多车场多车型车辆路径问题(Low-energy-consumption multi-depots heterogeneous-fleet vehicle routing problem with time windows, LMHFVPR_TW), 提出一种结合聚类分解策略的增强蚁群算法(Enhanced ant colony optimization based on clustering decomposition, EACO_CD)进行求解. 首先, 由于该问题具有强约束、大规模和NP-Hard等复杂性, 为有效控制问题的求解规模并合理引导算法在优质解区域搜索, 根据问题特点设计两种基于K-means的聚类策略, 将LMHFVPR_TW合理分解为一系列带时间窗的低能耗单车场单车型车辆路径子问题(Low-energy-consumption vehicle routing problem with time windows, LVRP_TW); 其次, 本文提出一种增强蚁群算法(Enhanced ant colony optimization, EACO)求解分解后的各子问题(LVRP_TW), 进而获得原问题的解. EACO不仅引入信息素挥发系数控制因子进一步动态调节信息素挥发系数, 从而有效控制信息素的挥发以提高算法的全局搜索能力, 而且设计基于4种变邻域操作的两阶段变邻域局部搜索(Two-stage variable neighborhood search, TVNS)来增强算法的局部搜索能力. 最后, 在不同规模问题上的仿真和对比实验验证了所提EACO_CD的有效性.Abstract: In this paper, an enhanced ant colony optimization combined with clustering decomposition strategy (EACO_CD) is proposed for solving the low-energy-consumption multi-depots heterogeneous-fleet vehicle routing problem with time windows (LMHFVPR_TW). Firstly, since the considered problem is a complex one with strong constraints, large scale and NP-hardness, In order to control the scale of problem and reasonably guide the algorithm to search in the high-quality solution region, two kinds of clustering methods based on K-means strategies are designed to reasonably decompose it into a series of subproblems (i.e., the low-energy-consumption vehicle routing problems with time windows (LVRP_TW) by utilizing the problem characteristics. Secondly, an enhanced ant colony optimization (EACO) to solve the decomposed subproblems (LVRP_TW) for obtaining the solution of the original problem is proposed. Not only a control factor of pheromone decay parameter in EACO is added to adjust the pheromone decay parameter dynamically so as to control the volatilization of pheromone effectively and improve the global search ability of ACO, but also a two-stage variable neighborhood search (TVNS) is designed based on four variable neighborhood operations to enhance its the local search ability. Finally, simulation experiments and comparisons on instances with different scales demonstrate the effectiveness of proposed EACO_CD.
-
近年来, 多智能体系统(Multi-agent systems, MASs)分布式协同控制问题的研究取得了显著进展, 引发各个领域的广泛关注. 该研究范畴涵盖生物系统中的群体行为[1]、分布式传感器网络技术[2]和智能电网管理[3]等多个方面. 一致性问题作为支撑MASs分布式协同控制的基础问题, 不仅在理论层面具有深远的意义, 而且在实际应用中展现出巨大的价值. 一致性控制的根本挑战在于设计高效的一致性算法或协议, 旨在确保MASs的所有智能体能够逐步调整其状态或输出, 最终达到相同, 即实现智能体的一致性.
目前, MASs一致性控制的研究可以根据系统中领航者的数量划分为三种类别: 无领航者的一致性控制、领导−跟随一致性控制(一个领航者)以及包含多个领航者的一致性控制[4−6]. 到目前为止, MASs一致性控制的研究涵盖越来越复杂的智能体动态特性和通信网络拓扑, 包括但不限于线性[7]或非线性MASs[8]、整数阶[9]或分数阶模型[10]、固定[11]或时变拓扑[12]、输入延迟[13]、输入饱和[14]等. 在上述复杂情况下, 各种适当的控制算法被提出以实现一致性控制. 此外, 由于智能体的通信和计算资源有限, 基于事件触发的控制策略[15−16]被用于实现一致性控制, 有效减少了不必要的能源消耗. 然而, 这些研究成果只能实现渐近一致性控制, 即在理论上调节时间趋于无穷. 在实际应用中, 由于渐近一致性的收敛时间较长, 难以满足任务的时效性需求.
相比之下, 有限时域一致性被认为是一种更为理想的控制策略. 有限时域控制不仅能够缩短闭环系统的收敛时间, 还具备更好的鲁棒性和抗干扰能力[17]. 文献[18]提出一种分散模型预测控制方案, 实现了一阶MASs的有限时域状态一致性控制. 文献[19]采用分布式线性二次型博弈方法, 实现了离散时间二阶MASs的有限时域状态一致性控制. 此外, 文献[20−23]研究离散时变MASs的$ H_\infty $有限时域状态一致性控制问题. 上述有限时域状态一致性协议的设计通常假设智能体动力学模型已知[20−23], 或仅考虑简单的一阶[18]、二阶[19]系统. 然而, 一阶和二阶系统无法充分描述实际系统的动态特性, 而且在实际应用中, 系统模型信息通常是未知的或难以获取的. 传统的有限时域一致性协议在系统模型未知的情况下并不适用, 难以满足实际应用的需求.
自适应动态规划(Adaptive dynamic programming, ADP)[24]或强化学习(Reinforcement learning, RL)[25]能够利用仿生学习机制解决系统模型未知情况下的优化控制问题[26]. 其中, 学习状态−动作值函数的Q学习算法[27]为实现无模型最优控制提供了一种有效的解决方案. 近年来, 学者们利用ADP或RL算法, 通过逼近求解耦合的哈密顿−雅可比−贝尔曼(Hamilton-Jacobi-Bellman, HJB)方程, 以实现MASs的最优渐近一致性控制[7−8, 28−31]. 例如, 基于Q学习的算法已经应用于异构MASs[7−8, 28−29]和同构MASs[30−31]中, 用以实现模型无关的最优一致性控制. 然而, 这些文献主要关注无限时域一致性控制问题. 相比之下, 模型无关的有限时域一致性控制问题更具挑战性, 因为它需要在满足值函数终端约束条件的同时求解耦合的时变HJB方程.
为解决上述问题, 学者们开始研究基于ADP或RL的算法, 以逼近耦合的时变HJB方程的近似解, 从而实现MASs的有限时域最优一致性控制. 文献[32]提出一种基于局部动力学的离策略(Off-policy) RL算法, 实现线性MASs的有限时域最优状态一致性控制. 此外, 文献[33]针对非线性MASs提出基于ADP的有限时域鲁棒事件触发最优状态一致性控制方法. 然而, 上述一致性控制器的设计[32−33]仍然依赖于MASs的部分模型信息, 而在实际情况下, 这些系统模型信息通常难以获得.
为克服系统模型必须已知的问题, 文献[34]采用神经网络逼近每个智能体的动态特性, 然后在神经网络模型的基础上基于ADP设计有限时域最优编队控制方法. 然而, 这种方式会产生额外的计算开销, 并引入逼近误差, 从而影响ADP方法的有效性. 文献[35]提出一种基于积分RL算法和零和博弈理论的模型无关有限时域鲁棒最优编队包含控制方法.
由于在实际的MASs中普遍存在执行器饱和的问题, 如无人车电机的输出转矩受最大功率限制, 无人机的舵面受物理结构限制等, 饱和的非线性特性通常会导致系统性能下降, 甚至可能导致系统不稳定, 使得执行器饱和问题在理论和实践上都极具挑战性. 上述研究结果[32−35]无法确保在模型未知的情况下实现具有执行器饱和约束的MASs一致性控制.
为解决这一问题, 学者们提出基于RL或ADP的方法来处理执行器饱和的MASs模型无关一致性控制问题. 例如, 文献[36]提出一种新型的辨识−评价−执行结构, 结合粘性消失法, 解决了有输入约束MASs的领导−跟随最优一致性控制问题. 文献[37]提出一种离策略RL算法, 通过逼近求解具有非二次型代价函数的耦合HJB方程, 以实现一致性控制. 文献[31, 38−40]使用低增益反馈(Low gain feedback, LGF)方法[41]处理执行器饱和问题, 并结合ADP方法实现执行器饱和的MASs最优一致性控制. 然而, 这些基于ADP或RL的模型无关一致性控制方法主要解决的是存在执行器饱和的MASs无限时域一致性控制问题, 只能实现渐近一致性控制, 即理论调控时间趋于无穷. 文献[42]基于ADP研究具有对称或不对称输入约束条件的MASs事件驱动有限时域最优状态一致性控制问题, 但其控制器的设计要求已知系统的模型信息.
受上述分析的启发, 本文将LGF方法与Q学习相结合, 用以解决执行器饱和的离散时间线性MASs模型无关有限时域一致性控制问题. 首先, 根据LGF方法的思想, 推导得到修正的时变黎卡提方程(Modified time-varying Riccati equation, MTVRE). 求解MTVRE可以得到时变的低增益反馈律, 同时可以通过调整低增益参数来避免执行器饱和. 然后, 参考文献[43−44], 设计依赖于系统状态、控制输入和低增益参数的时变参数化Q函数(Time-varying parameterized Q-function, TVPQF). 在TVPQF的基础上, 提出一种基于Q学习后向时间迭代模型无关一致性控制方法, 该方法在不需要已知系统动力学模型的前提下, 能够逼近求解MTVRE, 从而实现离散时间MASs的有限时域一致性控制.
本文将LGF方法与Q学习相结合, 提出一种针对执行器饱和的模型无关有限时域一致性控制方法. 主要贡献如下: 设计一种依赖于智能体状态、控制输入和低增益参数的TVPQF. 基于TVPQF, LGF控制器的设计减少了对系统动力学模型的依赖; 提出一种可以迭代更新低增益参数的后向时间模型无关控制算法, 并证明所提算法得到的时变LGF控制增益矩阵收敛于MTVRE的解; 另外, 证明所提算法不仅可以实现半全局一致性, 而且可以保证执行器饱和条件下的全局一致性控制, 并通过仿真实验进行论证.
本文的结构安排如下: 第1节首先介绍代数图论的相关知识, 并结合LGF方法介绍执行器饱和的离散时间MASs有限时域一致性控制问题的基于模型的求解方案. 第2节首先证明可以将执行器饱和的离散时间MASs有限时域一致性控制问题转化为执行器饱和的单智能体的最优控制问题, 接着提出基于TVPQF的后向时间迭代算法以逼近求解最优控制问题对应的MTVRE. 第3节提供仿真结果验证本文方法的有效性, 并进行对比实验, 比较性能指标突显本文方法的优越性. 第4节为结束语.
符号说明: $ {\bf{R}} $表示实数集, $ {\bf{R}}^{n \times m} $表示$ n\times m $维矩阵. $ I $表示具有兼容维数的单位矩阵. $ {\bf{0}} $表示具有兼容维数的全零向量或矩阵. $ \lambda_i(A) $表示矩阵$ A $的第$ i $个特征值. $ \text{Re} $表示实部. $ \text{rank}(A) $表示矩阵$ A $的秩. $ \text{argmax} $表示最大值索引. $ \text{argmin} $表示最小值索引. $ \text{vec} $为矩阵的拉直运算, 把矩阵按照列的顺序一列接一列的组成一个长向量. $ x^\text{T} $表示向量$ x $的转置.
1. 预备知识
1.1 代数图论
有$ N $个节点的加权图可记为$ G = (V,\; E,\; D) $, 其中节点和边的集合记为$ V = \{v_1,\;v_2,\;\cdots,\; v_N \} $和$ E = \{(v_i,\; v_j): v_i,\; v_j\in V\} $. 节点之间的连接关系由行随机矩阵$ D = [d_{ij}]\in {\bf{R}}^{N \times N} $决定, 其中$ d_{ii} > 0 $, $ \sum_{j = 1}^{N}d_{ij} = 1 $. 如果$ (v_i,\; v_j)\in E $, $ d_{ij}> 0 $; 否则$ d_{ij} = 0 $. 对于无向图$ G $, 行随机矩阵$ D $是对称的, 如果在任何一对不同的节点之间存在一条路径, 则称无向图$ G $是连通的. $ I-D $可看作是一种特殊的拉普拉斯矩阵, 满足$ \text{Re}(\lambda_1(I-D))< \text{Re}(\lambda_2(I-D))\le \cdots \le \text{Re}(\lambda_N(I-D)) $. 此外, 当且仅当有向图$ G $包含一个有向生成树, 或无向图$ G $连通时, 1是$ D $的一个单特征值. 令$ r\in {\bf{R}}^N $表示与$ I-D $的零特征值相关的左特征向量, 其满足$ r^\text{T} {\bf 1} = 1 $.
1.2 问题描述
考虑由$ N $个执行器饱和的智能体组成的离散时间MASs:
$$ \begin{equation} x_i(k+1) = Ax_i(k)+B\varrho (u_i(k)),\; i = 1,\;2,\;\cdots,\;N \end{equation} $$ (1) 式中, $ x_i(k) \in {\bf{R}}^n $, $ u_i(k) \in {\bf{R}}^m $分别表示智能体$ i $的状态向量以及输入向量; $ \varrho (\cdot):{\bf{R}}^m\rightarrow {\bf{R}}^m $表示饱和函数, 对于$ j = 1,\;2,\;\cdots,\;m $满足:
$$ \varrho (u_i^j(k)) = \left\{ \begin{aligned} & -c, & & \,u_i^j(k)<-c\\ & \, u_i^j(k), & &-c\leq u_i^j(k)\leq c\\ &\, c, & &\,u_i^j(k) > c \end{aligned} \right. $$ (2) 式中, $ c $表示饱和极限.
假设 1. 本文中, 智能体的系统模型是确定且未知的, 即$ A \in {\bf{R}}^{n \times n} $, $ B \in {\bf{R}}^{n \times m} $表示确定性的未知系统矩阵.
假设 2. 系统矩阵$ \left(A,\; B\right) $为输入有界下渐近零可控(Asymptotically null controllable with bounded controls, ANCBC), 即系统$ \left(A,\; B\right) $是可控的, 且$ A $的所有特征值都在单位圆上或单位圆内[41].
假设 3. 本文所考虑的用以描述离散时间MASs (1)拓扑结构的无向图$ G $是连通的.
假设 4. 本文所考虑的离散时间MASs (1)的阶次已知, 即$ n $是已知的.
本文研究的是具有执行器饱和的离散时间MASs的有限时域一致性控制问题. 所考虑的具体问题是: 在有限的时间区间内, 通过适当的控制策略设计, 使得所有智能体的状态在终端时刻达到一致, 即$ \lim _{k \rightarrow \tau}\left\|x_i(k)-x_j(k)\right\| = 0 $. 这种有限时域一致性控制要求在给定的时间范围$ \tau $内, 使所有智能体的状态在终端时刻达到某个共同的期望状态, 而不是在无限时域上渐近趋于一致.
参考文献[31], 针对离散时间MASs (1)可以设计如下状态反馈控制律:
$$ \begin{equation} u_i(k) = K(k) \sum\limits_{j = 1}^N d_{i j}\left(x_i\left(k\right)-x_j\left(k\right)\right) \end{equation} $$ (3) 其中, $ K(k) $为待设计的反馈控制增益矩阵.
引理 1. 对于具有$ N $个节点的离散时间MASs (1), 如果其对应的无向图$ G $是联通的, 则有$ \mu = 4 / (N(N - 1)) \leq \lambda_2 \left(I-D\right) $[45].
引理 2. 如果假设2和假设3成立, 则对于给定的有界集$ {\cal{X}} \in {\bf{R}}^n $, $ \forall x_i(0)\in {\cal{X}},\; i = 1,\;2,\; \cdots,\; N $, 存在最优低增益参数$ \varepsilon ^*\in (0,\;1] $, 对于任意$ \varepsilon \in (0, \varepsilon ^*] $, 离散时间MASs (1)可以在控制协议(3)下实现半全局一致性, 其中最优反馈控制增益矩阵满足:
$$ \begin{equation} K_{\varepsilon}^*(k) = -\left(B^{\mathrm{T}}P^*_{\varepsilon}(k+1)B+I\right)^{-1}B^{\mathrm{T}}P^*_{\varepsilon}(k+1)A \end{equation} $$ (4) 式中, $ P^*_{\varepsilon}(k) $满足如下MTVRE:
$$ \begin{split} P^*_{\varepsilon}(k) =\;& A^{\mathrm{T}} P^*_{\varepsilon}(k+1) A+\varepsilon I-\left(2 \mu-\mu^2\right) \;\times\\ & A^{\mathrm{T}}P^*_{\varepsilon}(k+1) B(B^{\mathrm{T}} P^*_{\varepsilon}(k+1) B+I)^{-1} \;\times \\ & B^{\mathrm{T}} P^*_{\varepsilon}(k+1) A\\[-1pt] \end{split} $$ (5) 同时, $ \lim _{\varepsilon \rightarrow 0}P^*_{\varepsilon}(k) = 0 $是单调的[38].
注1. 文献[38]考虑的是无限时域MASs一致性控制问题, 需求解修正的时变黎卡提方程. 而本文考虑的是有限时域一致性控制问题, 需求解MTVRE (5), 得到的正定矩阵$ P^*_{\varepsilon}(k) $以及LGF矩阵$ K_{\varepsilon}^*(k) $是时变的. 同时, 结合文献[38]中的引理2以及文献[46], 容易推导得到引理2.
注2. 相比于式(3), 式(4)中$ K_{\varepsilon}^*(k) $加下标$ \varepsilon $的原因在于, 根据LGF方法的思想, 反馈控制增益矩阵$ K_{\varepsilon}^*(k) $可以通过$ \varepsilon $进行调整, 从而使控制输入满足执行器饱和约束.
由引理2可知, 求解MTVRE (5)需要已知系统的模型参数$ \left(A,\; B\right) $. 然而, 在实际应用中, 系统的精确模型信息往往难以获取, 即便通过系统辨识可以获得模型信息, 但不可避免地会引入辨识误差. 同时, 引理2中给出的求解MTVRE (5)的方法只能实现半全局一致性. 为了解决这一问题, 本文首先将MASs的有限时域一致性控制问题转化为单智能体的有限时域最优控制问题, 并在无需系统模型信息且不依赖系统辨识的前提下, 提出一种结合低增益反馈与Q学习的模型无关数据驱动控制方法. 该方法能够在面对执行器饱和的情况下, 动态调整低增益参数, 从而在任意给定的智能体初始状态下, 实现离散时间MASs (1)的有限时域全局一致性控制.
2. 数据驱动有限时域一致性控制
在本节中, 将首先介绍使用LGF方法求解执行器饱和的单个智能体的优化控制问题, 进而推导得到MTVRE (5). 然后, 将介绍如何利用数据驱动方法, 通过单个智能体的可测量数据, 在系统模型信息未知的情况下, 逼近LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, 从而实现离散时间执行器饱和MASs (1)的有限时域一致性控制.
2.1 执行器饱和的单智能体优化控制方法
考虑如下执行器饱和的离散时间系统:
$$ \begin{equation} x_i(k+1) = Ax_i(k)+B\varrho (\zeta _i(k)) \end{equation} $$ (6) 其中, $ \zeta _i(k) $表示新的控制输入. 在接下来的基于Q学习的算法中, 将使用它来学习LGF矩阵 $ K_{\varepsilon}(k) $.
定义如下有限时域性能指标:
$$ \begin{equation} \begin{aligned} J_i = \sum\limits_{k = 0}^{\tau -1}r_{i}\left(x_{i}(k),\;\zeta _{i}(k),\;\varepsilon\right)+\varepsilon x_{i}^{\mathrm{T}}(\tau)x_{i}(\tau) \end{aligned} \end{equation} $$ (7) 式中, 最后一项$ \varepsilon x_{i}^{\mathrm{T}}(\tau)x_{i}(\tau) $代表终端约束条件; $ r_{i}\left(x_{i}(k),\;\zeta _{i}(k),\;\varepsilon\right) $表示智能体 $ i $的效用函数:
$$ \begin{equation} r_{i}\left(x_{i}(k),\;\zeta _{i}(k),\;\varepsilon\right) = \varepsilon x^{\mathrm{T}}_{i}(k)x_{i}(k)+\zeta^{\mathrm{T}}_{i}(k)\zeta_{i}(k) \end{equation} $$ (8) 根据有限时域性能指标 (7), 每个智能体$ i $的值函数可以表示为:
$$ \begin{equation} V_i(x_{i}(k)) = \sum\limits_{j = k}^{\tau-1} r_{i}\left(x_{i}(j),\;\zeta _{i}(j),\;\varepsilon\right)+\varepsilon x_{i}^{\mathrm{T}}(\tau)x_{i}(\tau) \end{equation} $$ (9) 下面引理证明, 当控制输入$ \zeta _i(k) = \mu K_{\varepsilon}(k)\;\times x_i(k) $时, 值函数(9)可以表示为二次型形式.
引理3. 如果离散时间系统(6)的控制输入可以表示为$ \zeta _i(k) = \mu K_{\varepsilon}(k)x_i(k) $, 则智能体$ i $的值函数$ V_i(x_{i}(k)) $可以表示为如下二次型形式:
$$ \begin{equation} V_i(x_{i}(k)) = x_{i}^{\mathrm{T}}(k) P_\varepsilon(k) x_{i}(k) \end{equation} $$ (10) 式中, $ P_\varepsilon(k) = P^{\mathrm{T}}_\varepsilon(k)>0 $. $ P_\varepsilon(\tau) = \varepsilon I $.
证明. 本部分将基于最优性原理, 利用终端约束条件采用后向时间的方式进行证明.
当$ k = \tau $时, 可以很容易地从式(9)得到:
$$ \begin{equation} V_i(x_{i}(\tau)) = \varepsilon x_{i}^{\mathrm{T}}(\tau)x_{i}(\tau) \end{equation} $$ (11) 因此, 可以得到$ P_\varepsilon (\tau) = P^{\mathrm{T}}_\varepsilon (\tau) = \varepsilon I $.
当$ k = \tau -1 $时, 结合式(8)和(9)可以得到:
$$ \begin{split} V_i(x_{i}(\tau -1)) =\; &\varepsilon x^{\mathrm{T}}_{i}(\tau -1)x_{i}(\tau -1)\;+\\ &\zeta^{\mathrm{T}}_{i}(\tau -1)\zeta_{i}(\tau -1) +\varepsilon x_{i}^{\mathrm{T}}(\tau)x_{i}(\tau) \end{split} $$ (12) 将式(6)代入式(12)中, 得到:
$$ \begin{split} V_i(x_{i}(\tau -1)) =\; &\varepsilon x^{\mathrm{T}}_{i}(\tau -1)x_{i}(\tau -1)\;+\\ &\zeta^{\mathrm{T}}_{i}(\tau -1)\zeta_{i}(\tau -1)\;+\\ &\varepsilon(Ax_i(\tau -1)+B\zeta _i(\tau -1))^{\mathrm{T}}\;\times \\ &(Ax_i(\tau -1)+B\zeta _i(\tau -1))\\[-1pt] \end{split} $$ (13) 然后, 将$ \zeta _i(\tau -1) = \mu K_{\varepsilon}(\tau -1)x_i(\tau -1) $代入式(13)中, 可以得到:
$$ \begin{split} V_i(x_{i}(\tau -1)) = \; &x^{\mathrm{T}}_{i}(\tau -1)[\varepsilon I + \mu^2 K_{\varepsilon}^{\mathrm{T}}(\tau -1)K_{\varepsilon}(\tau \;-\\ &1) +\varepsilon\left(A + \mu BK_{\varepsilon}(\tau -1)\right)^{\mathrm{T}}\;\times\\ &(A + \mu BK_{\varepsilon}(\tau -1))]x_{i}(\tau -1) \\[-1pt]\end{split} $$ (14) 当$ k = \tau -1 $时, 从式(14)可以得到:
$$ \begin{split} P_\varepsilon (\tau -1) = \; &\varepsilon I + \mu^2 K_{\varepsilon}^{\mathrm{T}}(\tau -1)K_{\varepsilon}(\tau -1)\;+ \\ &\varepsilon\left(A + \mu BK_{\varepsilon}(\tau -1)\right)^{\mathrm{T}}\;\times\\ &(A + \mu BK_{\varepsilon}(\tau -1)) \end{split} $$ (15) 从上式可以得到 $ P_\varepsilon (\tau -1) = P^{\mathrm{T}}_\varepsilon (\tau -1)>0 $.
采用与$ P_\varepsilon (\tau -1) $相同的方式, 可以类似地确定, 对于$ k = 0,\;1,\; \cdots,\; \tau-2 $, 矩阵$ P_\varepsilon(k) $也符合$ P_\varepsilon(k) = P^{\mathrm{T}}_\varepsilon (k)>0 $.
□ 下面定理将证明, 针对执行器饱和的离散时间系统(6)以及对应的有限时域性能指标(7), 存在最优的LGF控制增益矩阵$ K_{\varepsilon}(k) $, 使得智能体$ i $的值函数$ V_i(x_{i}(k)) $可以表示为式(10).
定理1. 考虑执行器饱和离散时间系统(6)以及对应的有限时域性能指标(7), 其最优控制律满足:
$$ \begin{equation} {\zeta}^*_i(k) = K^*_{\varepsilon}(k)x_i(k) \end{equation} $$ (16) 其中, $ K^*_{\varepsilon}(k) $满足式(4). 如果令$ {\zeta}^*_i(k) = \mu K^*_{\varepsilon}(k)x_i(k) $, 则$ P^*_{\varepsilon}(k) $满足式(5).
证明. 根据值函数的定义(9)可知, 值函数满足如下贝尔曼方程:
$$ \begin{equation} \begin{aligned} V_i(x_{i}(k)) = \varepsilon x^{\mathrm{T}}_{i}(k)x_{i}(k)+{\zeta}^{\mathrm{T}}_{i}(k){\zeta}_{i}(k)+V_i(x_{i}(k+1)) \end{aligned} \end{equation} $$ (17) 同时, 最优值函数满足:
$$ \begin{split} V_i^*(x_i(k)) =\;&\min_{{\zeta}_i(k)} \sum\limits_{j = k}^{\tau-1}(\varepsilon x^{\mathrm{T}}_{i}(j)x_{i}(j)+{\zeta}^{\mathrm{T}}_{i}(j){\zeta}_{i}(j)\;+\\ &\varepsilon x_{i}^{\mathrm{T}}(\tau )x_{i}(\tau)) \\[-1pt]\end{split} $$ (18) 结合式(17)和式(18), 可以得到如下的贝尔曼最优方程:
$$ \begin{split} V^*_i(x_{i}(k)) =\; &\min_{{\zeta}_i(k)}\left(\varepsilon x_{i}^{\mathrm{T}}(k)x_{i}(k)+{\zeta}_i^{\mathrm{T}}(k){\zeta}_i(k)\right.+\\ &\left.V^*_i(x_{i}(k+1))\right) \end{split} $$ (19) 当$ k = \tau-1 $时, 由式(18)可知:
$$ \begin{split} V^*_i(x_{i}(\tau-1)) = \;&\min_{{\zeta}_i(\tau-1)} \left(\varepsilon x^{\mathrm{T}}_{i}(\tau-1)x_{i}(\tau-1)\right.+\\ &\left.{\zeta}^{\mathrm{T}}_{i}(\tau-1){\zeta}_{i}(\tau-1)+\varepsilon x_{i}^{\mathrm{T}}(\tau )x_{i}(\tau)\right) \end{split} $$ (20) 将式(6)代入式(20)中, 得到:
$$ \begin{split} V^*_i(x_{i}(\tau-1)) = \; &\min_{{\zeta}_i(\tau-1)}(\varepsilon x_{i}^{\mathrm{T}}(\tau-1)x_{i}(\tau-1)\;+\\ &{\zeta}_i^{\mathrm{T}}(\tau-1){\zeta}_i(\tau-1)\;+\\ &\varepsilon(Ax_{i}(\tau-1)+B{\zeta}_i(\tau-1))^{\mathrm{T}}\;\times\\ &(Ax_{i}(\tau-1)+B{\zeta}_i(\tau-1))) \\[-1pt] \end{split} $$ (21) 从式(21)可以得到最优控制策略满足:
$$ \begin{split} {\zeta}^*_i(\tau-1) = \;&\arg\min_{{\zeta}_i(\tau-1)}(\varepsilon x_{i}^{\mathrm{T}}(\tau-1)x_{i}(\tau-1)\;+\\ &{\zeta}_i^{\mathrm{T}}(\tau-1){\zeta}_i(\tau-1)\;+\\ &\varepsilon(Ax_{i}(\tau-1)+B{\zeta}_i(\tau-1))^{\mathrm{T}}\;\times\\ & (Ax_{i}(\tau-1)+B{\zeta}_i(\tau-1))) \end{split} $$ (22) 为了得到最优控制策略, 可以通过上式右半部分对$ {\zeta}_i(\tau-1) $求导, 并令导数为零. 则有:
$$ \begin{equation} 2{\zeta}^{\mathrm{T}}_i(\tau-1)+2\varepsilon\left(Ax_{i}(\tau-1)+B{\zeta}_i(\tau-1)\right)^{\mathrm{T}}B = 0 \end{equation} $$ (23) 因此, 可以得到最优控制策略:
$$ {\zeta}^*_i(\tau-1) = -\varepsilon\left(\varepsilon B^{\mathrm{T}}B+I\right)^{-1}B^{\mathrm{T}}Ax_i(\tau-1) $$ (24) 结合值函数的终端约束条件可知$ P^*_{\varepsilon}(\tau) = \varepsilon I $, 则式(24)可以重写为:
$$ \begin{split} {\zeta}^*_i(\tau-1) = \; &\left(B^{\mathrm{T}} P^*_{\varepsilon}(\tau) B+I\right)^{-1}\times\\ &B^{\mathrm{T}}P^*_{\varepsilon}(\tau) Ax_i(\tau-1) = \\ &K^*_{\varepsilon}(\tau-1)x_i(\tau-1) \end{split} $$ (25) 比较式(4)和式(25), 可知$ K^*_{\varepsilon}(\tau-1) $满足式(4).
结合文献[46]以及引理3, 可以得到最优值函数$ V^*_i(x_{i}(\tau-1)) $可以写成如下形式:
$$ \begin{equation} V^*_i(x_{i}(\tau-1)) = x_{i}^{\mathrm{T}}(\tau-1) P^*_\varepsilon(\tau-1) x_{i}(\tau-1) \end{equation} $$ (26) 同时, 将$ {\zeta}^*_i(\tau-1) = \mu K^*_{\varepsilon}(\tau-1)x_i(\tau-1) $代入式(20)中, 很容易得到$ P^*_{\varepsilon}(\tau-1) $满足式(15).
采用与$ k = \tau -1 $相同的方式, 可以依次得到$ K^*_{\varepsilon}(k),\; k = \tau -2,\; \cdots,\; 1,\; 0 $满足式(4), 并且值函数满足:
$$ \begin{equation} V^*_i(x_{i}(k)) = x_{i}^{\mathrm{T}}(k) P^*_\varepsilon(k) x_{i}(k) \end{equation} $$ (27) 此外, 将$ {\zeta}^*_i(k) = \mu K^*_{\varepsilon}(k)x_i(k) $代入式(19)中, 并结合式(27), 很容易得到$ P^*_{\varepsilon}(k) $满足式(15).
□ 注3. 与有限/固定时间控制不同, 本文所考虑的有限时域一致性控制是指控制器在一个预算的时间段内进行设计. 在这个时间段结束时, 控制器的目标是使系统状态达到某个期望的状态或者满足特定的性能指标. 有限时域控制问题通常涉及优化一个性能指标函数, 该函数定义在从初始时刻到终止时刻的时间区间上, 如本文所考虑的有限时域性能指标函数(7), 并且需要考虑在此期间系统的动态行为和可能存在的约束条件, 如本文所考虑的执行器饱和约束. 而有限时间控制强调的是收敛时间$ t $趋于一个固定值$ T $达到稳定, 该$ T $是根据初值和控制参数计算出来的. 固定时间控制是一种特殊的有限时间控制, 也是$ t $趋于一个固定值$ T $达到稳定, 该$ T $的计算和初值无关, 但是计算的$ T $有保守性. 有限时域控制可以看作有限时间控制的一种特殊情况, 其侧重点在于需要在固定时间范围内优化一个性能指标函数.
注4. 根据低增益反馈控制方法[41], 可以对低增益参数进行动态调整, 逐步将控制输入限制在饱和值范围内, 从而避免执行器饱和现象. 在引理3以及定理1的证明过程中, 由于低增益参数的存在, 在证明过程中假定通过调整低增益参数得到满足执行器饱和约束的控制输入. 因此, 在涉及控制输入的证明过程中, 饱和函数$ \rho ( \cdot ) $没有显示地出现.
从以上分析可知, 可以将针对执行器饱和的离散时间MASs (1)的有限时域一致性控制问题转化为针对执行器饱和的离散时间系统(6)以及有限时域性能指标(7)的最优控制问题. 不同之处在于, 为了实现有限时域一致性控制, 需要改变由最优控制问题求得的控制策略. 同时, 依据LGF方法的特点, 可以通过调整低增益参数$ \varepsilon $实现避免执行器饱和的目标.
2.2 模型无关有限时域控制方法
在这一部分, 首先, 结合Q学习的思想定义TVPQF; 然后, 提出一种数据驱动的后向时间迭代方法, 在仅需要单个智能体可测量数据的前提下, 逼近求解MTVRE (5), 以实现有限时域一致性控制.
依据文献[27], 定义如下TVPQF:
$$ \begin{split} &Q_{\varepsilon}\left(x_{i}(k),\;\zeta _{i}(k),\;\tau-k\right) = \\ &\qquad r_{i}\left(x_{i}(k),\;\zeta _{i}(k),\;\varepsilon\right)+V_{i}^*\left(x_{i}(k+1)\right) \end{split} $$ (28) 其中, $ Q_{\varepsilon}(x_{i}(\tau)) = \varepsilon x^{\mathrm{T}}_{i }(\tau)x_{i}(\tau) $.
定义变量$ \xi_{i}(k) = \left[x_{i}^{\mathrm{T}}(k),\; \zeta_{i}^{\mathrm{T}}(k)\right]^{\mathrm{T}} $. 同时, 将式(6)和(19)代入式(28)中, 可以得到:
$$ \begin{equation} \begin{aligned} Q_{\varepsilon}\left(x_{i}(k),\;\zeta _{i}(k),\;\tau-k\right) = \xi^{\mathrm{T}}_{i}(k){\cal{H}}_{\varepsilon}(k)\xi_{i}(k) \end{aligned} \end{equation} $$ (29) 式中, $ {\cal{H}}_{\varepsilon}(k) $表示TVPQF的核函数, 定义如下:
$$ \begin{split} &{{\cal{H}}_\varepsilon }(k): = \left[ {\begin{array}{*{20}{l}} {{{\cal{H}}_{xx}}(k)}&{{{\cal{H}}_{x\zeta }}(k)}\\ {{{\cal{H}}_{\zeta x}}(k)}&{{{\cal{H}}_{\zeta \zeta }}(k)} \end{array}} \right] = \\ &\;\;\;\;\left[ {\begin{array}{*{20}{c}} {\varepsilon I + {A^{\rm{T}}}{P_\varepsilon }(k + 1)A}&{{A^{\rm{T}}}{P_\varepsilon }(k + 1)B}\\ {{B^{\rm{T}}}{P_\varepsilon }(k + 1)A}&{{B^{\mathrm{T}}}{P_\varepsilon }(k + 1)B + I} \end{array}} \right] \end{split} $$ (30) 同时, 通过TVPQF的定义 (28)可以得到最优值函数与最优TVPQF的关系如下:
$$ \begin{split} V^*_i(x_i(k))=\; & \min_{{\zeta}_i(k)}Q_{\varepsilon}\left(x_{i}(k),\;\zeta _{i}(k),\;\tau-k\right) = \\ &Q^*_{\varepsilon}\left(x_{i}(k),\;\zeta^* _{i}(k),\;\tau-k\right) \end{split} $$ (31) 根据TVPQF的定义可知, 最优LGF控制律满足:
$$ \begin{equation} \zeta^*_i(k) = \arg\min\limits_{\zeta_i(k)}Q_{\varepsilon}\left(x_{i}(k),\;\zeta _{i}(k),\;\tau-k\right) \end{equation} $$ (32) 求解$ \frac{{\partial Q_{\varepsilon}(x_{i}(k),\;\zeta _{i}(k),\;\tau-k)} }{ {\partial \zeta_i(k)}} = 0 $, 可以得到:
$$ \begin{equation} K_{\varepsilon}^*(k) = -{\cal{H}}_{\zeta \zeta}^{*}{}^{-1}(k){\cal{H}}_{\zeta x}^*(k) \end{equation} $$ (33) 另外, 将$ {\zeta}^*_i(k) = \mu {K}^*_{\varepsilon}(k)x_i(k) $、式(33)代入式(31), 同时结合式(29), 得到:
$$ \begin{split} P^*_{\varepsilon}(k) =\; &{\cal{H}}^*_{x x}(k)-\mu{K}^*_{\varepsilon}(k){\cal{H}}^*_{\zeta x}(k)+\mu{\cal{H}}^*_{x \zeta}(k)\;\times\\ & {K}^{*,\;\mathrm{T}}_{\varepsilon}(k)+\mu^2{K}^*_{\varepsilon}(k){\cal{H}}^*_{\zeta \zeta}(k){K}^{*,\;\mathrm{T}}_{\varepsilon}(k)\\[-1pt] \end{split} $$ (34) 根据式(33)和(34)可知, 通过设计的TVPQF, 可以将计算$ P^*_{\varepsilon}(k) $转变为计算$ {\cal{H}}_{\varepsilon}^*(k) $, 以获取最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, 并且避免对系统模型信息的依赖. 下面将介绍如何采用后向时间的方式逼近求解$ {\cal{H}}_{\varepsilon}^*(k) $.
假设通过$ \eta $次实验, 收集到$ \eta $组样本数据$ \{x^j_{i}(k), \zeta^j _{i}(k),\;x^j_{i}(k+1)\} $, 其中$ j = 1,\;2,\;\cdots,\;\eta $.
当$ k = \tau -1 $时, 定义:
$$ \begin{split} {\cal{Q}}^j_{\varepsilon}(\tau - 1) = \;& \varepsilon x^{j,\;\mathrm{T}}_{i}(\tau - 1)x^j_{i}(\tau - 1)+\zeta^{j,\;\mathrm{T}}_{i}(\tau - 1)\;\times\\ & \zeta^j_{i}(\tau - 1)+x^{j,\; \mathrm{T}}_{i }(\tau)P^*_{\varepsilon}(\tau) x^j_{i}(\tau) \\[-1pt] \end{split} $$ (35) 式中, $ P^*_{\varepsilon}(\tau) = \varepsilon I $.
同时, 根据式(29), 可以得到$ {\cal{Q}}^j_{\varepsilon}(\tau - 1) $的另一种表达形式如下:
$$ \begin{equation} \hat{{\cal{Q}}}^j_{\varepsilon}(\tau - 1) = \xi^{j,\; \mathrm{T}}_{i}(\tau - 1){\cal{H}}_{\varepsilon}(\tau - 1)\xi^j_{i}(\tau - 1) \end{equation} $$ (36) 应用线性参数化方法, 式(36)可以重写成:
$$ \begin{equation} \hat{{\cal{Q}}}^j_{\varepsilon}(\tau - 1) = \bar{\xi}^{j,\; \mathrm{T}}_{i}(\tau - 1)\text{vec}({\cal{H}}_{\varepsilon}(\tau - 1)) \end{equation} $$ (37) 其中,
$$ \begin{split} \bar{\xi}^{j}_{i}(\tau - 1) =\;& [({\xi}^{1,\;j}_{i})^2,\; 2{\xi}^{1,\;j}_{i}{\xi}^{2,\;j}_{i},\; \cdots,\; 2{\xi}^{1,\;j}_{i}{\xi}^{l,\;j}_{i},\\ & ({\xi}^{2,\;j}_{i})^2,\;2{\xi}^{2,\;j}_{i}{\xi}^{3,\;j}_{i},\; \cdots,\;\\ &2{\xi}^{2,\;j}_{i}{\xi}^{l,\;j}_{i},\; \cdots,\; ({\xi}^{l,\;j}_{i})^2]^{\mathrm{T}} \nonumber \end{split} $$ 上面变量的表达式中$ l = n+m $表示变量$ \bar{\xi}^{j}_{i}(\tau \;- 1) $的维数. 另外, 为方便, 省去了$ \tau - 1 $.
结合式(35)和(37)可知, 可以通过求解如下优化方程用以获取TVPQF对应的最优核矩阵$ {\cal{H}}_{\varepsilon}^*(\tau - 1) $:
$$ \begin{split} \text{vec}({\cal{H}}_{\varepsilon}^*(\tau - 1)) = \; & \arg\min\sum\limits_{j = 1}^\eta(\bar{\xi}^{j,\; \mathrm{T}}_{i}(\tau - 1)\;\times\\ &\text{vec}({\cal{H}}_{\varepsilon}(\tau - 1))-{\cal{Q}}^j_{\varepsilon}(\tau - 1))^2 \end{split} $$ (38) 得到$ {\cal{H}}_{\varepsilon}^*(\tau - 1) $, 就可以通过式(33)求解最优LGF控制增益矩阵$ K_{\varepsilon}^*(\tau - 1) $, 以及通过式(34)获取最优值函数对应的核矩阵$ P_{\varepsilon}^*(\tau - 1) $.
依据求解$ {\cal{H}}_{\varepsilon}^*(\tau - 1) $的思路, 可以通过后向时间求解的方式逼近求解$ {\cal{H}}_{\varepsilon}^*(k) $, $ K_{\varepsilon}^*(k) $, 以及$ P_{\varepsilon}^*(k) $, $ k = \tau-2,\; \cdots,\; 1,\; 0 $.
当$ k = \tau-2,\; \cdots,\; 1,\; 0 $时, 定义:
$$ \begin{split} {\cal{Q}}^j_{\varepsilon}(k) =\; & \varepsilon x^{j,\;\mathrm{T}}_{i}(k)x^j_{i}(k)+\zeta^{j,\;\mathrm{T}}_{i}(k)\zeta^j_{i}(k)\;+\\ &x^{j,\; \mathrm{T}}_{i }(k+1)P^*_{\varepsilon}(k+1) x^j_{i}(k+1) \end{split} $$ (39) 同样地, 可以得到$ {\cal{Q}}^j_{\varepsilon}(k) $的另一种表达形式:
$$ \begin{equation} \hat{{\cal{Q}}}^j_{\varepsilon}(k) = \xi^{j,\; \mathrm{T}}_{i}(k){\cal{H}}_{\varepsilon}(k)\xi^j_{i}(k) \end{equation} $$ (40) 参照式(38), 可以得到如下优化问题:
$$ \begin{split} &\text{vec}({\cal{H}}_{\varepsilon}^*(k)) = \\ &\qquad \arg\min\sum\limits_{j = 1}^\eta\left(\bar{\xi}^{j,\; \mathrm{T}}_{i}(k)\text{vec}({\cal{H}}_{\varepsilon}(k))-{\cal{Q}}^j_{\varepsilon}(k)\right)^2 \end{split} $$ (41) 通过式(41)求解得到$ {\cal{H}}_{\varepsilon}^*(k) $, 就可以通过式(33)求解最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, 以及通过式(34)获取最优值函数对应的核矩阵$ P_{\varepsilon}^*(k) $. 下面将介绍如何求解优化问题(38)和问题(41). 由于两者具有相似性, 下面将问题(38)和问题(41)归结为一类问题进行介绍.
优化问题(38)和问题(41)可以写成如下形式:
$$ \begin{split} &\text{vec}\left({\cal{H}}_{\varepsilon}^*(k)\right) = \\ & \qquad\arg\min\sum\left(\bar{{\xi}}^{\,\mathrm{T}}_{i}(k)\text{vec}({\cal{H}}_{\varepsilon}(k))-{\cal{Q}}_{\varepsilon}(k)\right)^2 \end{split} $$ (42) 式中, $ \bar{{\xi}}_{i}(k) = \left[\bar{{\xi}}^{1}_{i}(k),\;\bar{{\xi}}^{2}_{i}(k),\;\cdots,\;\bar{{\xi}}^{\eta}_{i}(k)\right]^{\mathrm{T}} $; $ {\cal{Q}}_{\varepsilon}(k) = \left[{\cal{Q}}_{\varepsilon}^{1}(k),\; {\cal{Q}}_{\varepsilon}^{2}(k),\; \cdots,\; {\cal{Q}}_{\varepsilon}^{\eta}(k)\right]^{\mathrm{T}} $, $ k = 0,\;1,\;\cdots,\;\tau-1 $.
应用最小二乘法, 可以得到优化问题 (42)的解如下:
$$ \begin{equation} \text{vec}({\cal{H}}_{\varepsilon}^*(k)) = \left(\bar{{\xi}}_{i}(k)\bar{{\xi}}^{\,\mathrm{T}}_{i}(k)\right)^{-1}\bar{{\xi}}_{i}(k){\cal{Q}}_{\varepsilon}(k) \end{equation} $$ (43) 为确保优化问题(42)的解(43)的唯一性, 需要满足如下条件:
$$ \begin{equation} \text{rank}(\bar{{\xi}}_{i}(k)) = \frac{l(l+1)}{2} \end{equation} $$ (44) 即矩阵$ \bar{{\xi}}_{i}(k) $满秩.
如果搜集到的样本$ \{x^j_{i}(k),\;\zeta^j _{i}(k),\;x^j_{i}(k+1)\} $的数量$ \eta\ge {{l(l + 1)} / 2} $, 且每次实验收集到的数据之间服从高斯分布, 那么条件(44)成立[46].
由以上分析可知, 采用后向时间求解的方式可以得到最优TVPQF对应的核矩阵$ {\cal{H}}_{\varepsilon}^*(k) $. 同时, 由式(35)和(39)可知, TVPQF会受到低增益参数$ \varepsilon $的影响. 因此, 可以通过调整低增益参数$ \varepsilon $用以更新LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, 从而使控制器$ u_i(k) $避免输入饱和. 算法1对上面的论述进行了总结.
算法 1. 执行器饱和约束下模型无关有限时域一致性控制
输入. 实验次数$ \eta $, 低增益参数$ \varepsilon $, 有限时域$ \tau $.
输出. 最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, 以及最优TVPQF对应的核矩阵$ {\cal{H}}_{\varepsilon}^*(k) $, $ k = 0,\;1,\; \cdots,\; \tau -1 $.
1) 数据收集: 生成符合高斯分布的随机控制输入$ \{\zeta^j_{i}(0),\; \zeta^j_{i}(1),\; \cdots,\; \zeta^j_{i}(\tau-1)\} $, 以及随机初始状态变量$ x^j_{i}(0) $, $ j = 1,\;2,\;\cdots,\;\eta $, 应用于系统(6), 从而收集产生的样本数据$ \{x^j_{i}(k),\;\zeta^j _{i}(k),\;x^j_{i}(k+1)\} $, 其中$ j = 1,\;2,\;\cdots,\;\eta $; $ k = 0,\;1,\;\cdots,\;\tau-1 $.
2) 计算$ K_{\varepsilon}^*(\tau -1) $: 通过式(43)求解优化问题(38), 得到$ {\cal{H}}_{\varepsilon}^*(\tau -1) $. 结合式(33)推导得到最优LGF控制增益矩阵$ K_{\varepsilon}^*(\tau -1) $, 并将其存储.
3) 计算$ K_{\varepsilon}^*(k) $: 从$ k = \tau -2 $到$ k = 0 $, 依次通过式(43)求解优化问题(42), 迭代计算$ {\cal{H}}_{\varepsilon}^*(k) $. 结合式(33)推导得到最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, $ k = \tau\;- 2,\; \cdots,\; 1,\; 0 $, 并将其存储.
4) 饱和度检查: 对于每一个$ k = 0,\;1,\;\cdots,\;\tau-1 $, 验证
$$ \qquad\left\lVert u_i(k) \right\rVert_\infty = \left\lVert K_{\varepsilon}^*(k) \sum\limits_{j = 1}^N d_{i j}\left(x_i(k)-x_j(k)\right)\right\rVert_\infty \le c $$ 其中, $ i = 1,\;2,\;\cdots,\;N $. 如果不满足, 则减小$ \varepsilon $并重复步骤2)和步骤3).
5) 停止迭代: 当控制输入不再饱和时, 停止迭代过程.
注5. 算法1中, 低增益参数$ \varepsilon $可以通过比例规则进行调整: $ \varepsilon_{j+1} = \alpha \varepsilon_j $, 其中$ 0<\alpha<1 $. 另外需要强调的是, 控制输入的饱和度评估发生在其应用到MASs (1)之前. 因此, MASs在实际执行的过程中不会超过其执行器饱和约束.
注6. 算法1中的饱和度检查环节必然会受到智能体初始状态的影响, 不同的初始状态可能会最终得到不同的低增益参数$ \varepsilon $. 另外, 算法1的目的并不是寻找最优低增益参数$ \varepsilon^* $, 而是对于不同的初始状态寻找$ \varepsilon \in (0,\;\varepsilon^*] $, 从而得到对应的最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, 达到避免执行器饱和的目标.
下面定理将证明通过算法1得到的最优LGF控制增益矩阵是最优的.
定理2. 如果进行收集样本数据的实验次数$ \eta \ge {{l(l + 1)} / 2} $, 且收集得到的样本数据$ \{x^j_{i}(k),\, \zeta^j _{i}(k), x^j_{i}(k+1)\} $服从高斯分布, 则算法1得到的LGF控制增益矩阵$ K_{\varepsilon}^*(k) $是最优的, 也就是MTVRE (5)对应的解.
证明. 根据LGF方法的思想, 针对执行器饱和约束问题, 存在最优低增益参数$ \varepsilon^* $ [42]. 同时, 注意到算法1中关于低增益参数$ \varepsilon $的调整处于估计$ {\cal{H}}_{\varepsilon}(k) $的外循环. 因此, 低增益参数$ \varepsilon $不会影响TVPQF核矩阵$ {\cal{H}}_{\varepsilon}(k) $的收敛性. 假设低增益参数$ \varepsilon $在算法1中是固定的, 即考虑MTVRE (5)和TVPQF (23)中包含相同的低增益参数$ \varepsilon $的情况.
当初始样本数据$ x^j_{i}(0),\; j = 1,\;2,\;\cdots,\;\eta $, 以及$ \zeta^j _{i}(k),\; k = 0,\;1,\;\cdots,\;\tau-1 $服从高斯分布时, 很容易得到每次收集样本数据的实验是线性独立的. 此外, 如果实验次数$ \eta\ge {{l(l + 1)} / 2} $, 则式(43)中构造得到的数据矩阵$ \bar{{\xi}}_{i}(k) $, $ k = 0,\;1,\;\cdots,\;\tau-1 $满秩. 需要注意的是, 结合$ \text{vec}({\cal{H}}_{\varepsilon}^*(k)) $的定义以及式(30)可知,$ \text{vec}({\cal{H}}_{\varepsilon}^*(k)) $, $ k = 0,\;1,\;\cdots,\;\tau-1 $拥有$ {{l(l + 1)} / 2} $个独立元素. 结合矩阵$ \bar{{\xi}}_{i}(k) $满秩的结论, 可知优化问题(42)有唯一解, 即为式(43). 值得注意的是, 所设计的算法1以离线后向时间迭代的方式运行, 即利用终端约束条件$ P_{\varepsilon}(\tau) $从$ k = \tau -1 $开始依次向后计算$ \text{vec}({\cal{H}}_{\varepsilon}^*(k)) $. 同时, 式(43)构成了优化问题(42)的唯一解. 可以得出结论: 通过执行算法1得到的$ \text{vec}({\cal{H}}_{\varepsilon}^*(k)) $是最优的.
值得注意的是, 矩阵$ {\cal{H}}_{\varepsilon}^*(k) $是由$ {{l(l + 1)} / 2} $个元素组成的对称矩阵, $ \text{vec}({\cal{H}}_{\varepsilon}^*(k)) $表示矩阵$ {\cal{H}}_{\varepsilon}^*(k) $经过列排列之后组成的长向量. 由于算法1得到的$ \text{vec}({\cal{H}}_{\varepsilon}^*(k)) $是最优的. 因此, 算法1得到的结果$ {\cal{H}}_{\varepsilon}^*(k) $即为所定义的TVPQF的最优核矩阵. 结合式(33)以及引理2可知, 算法1得出的LGF控制增益矩阵$ K_{\varepsilon}^*(k) $也是最优的. 同时, 结合定理1以及式(31)可知, 通过算法1得到LGF控制增益矩阵$ K_{\varepsilon}^*(k) $等价于求解MTVRE (5).
□ 下面定理将证明算法1可以实现离散时间MASs (1)的全局有限时域一致性控制而不仅仅是半全局有限时域一致性控制.
定理3. 如果假设2和假设3成立, 通过算法1得到的LGF控制增益矩阵$ K_{\varepsilon}^*(k) $, 离散时间MASs (1)可以实现全局有限时域一致性控制.
证明. 如果假设2和假设3成立, 由引理2以及定理2可知, 算法1得到的LGF控制增益矩阵$ K_{\varepsilon}^*(k) $可以实现半全局有限时域一致性控制. 在算法1中, 如果控制输入违反执行器饱和, 则在下次迭代时会减小低增益参数$ \varepsilon $, 因此必然可以找到一个足够小的$ \varepsilon \in (0,\;\varepsilon^*] $满足执行器饱和. 另外, 从定理2的证明过程可知, 如果$ \varepsilon $固定, 由算法1得到的TVPQF核矩阵以及LGF控制增益矩阵均是最优的, 且可以实现有限时域最优一致性控制. 如果智能体的初始状态不同, 必然会迭代得到一个固定的低增益参数$ \varepsilon $, 对应地, 即可通过算法1得到LGF控制增益矩阵. 因此, 算法1可以实现离散时间MASs (1)的全局有限时域一致性控制.
□ 3. 仿真实验
本节首先建立一个仿真实验, 来说明本文方法的有效性; 然后进行对比实验, 用本文方法与对比方法进行仿真实验, 用评价指标结果说明本文方法的优越性.
3.1 仿真实验1
考虑一个由6个智能体组成的离散时间MASs, 其动力学方程为(1), 相关的矩阵为:
$$ \begin{equation} A = \begin{bmatrix} 0 & 1\\ -1 & 1 \\ \end{bmatrix},\;\quad B = \begin{bmatrix} -1 \\ 0 \end{bmatrix} \end{equation} $$ (45) 矩阵$ A $的特征值$ 0.5\pm0.866\mathrm{i} $都在单位圆内, 且$ (A,\; B) $是可控的. 因此, 假设2成立. 在本节仿真中, 执行器饱和函数的饱和阈值设为$ c = 1 $. 离散时间MASs的通信拓扑可以用图1所示的无向图表示. 从图中可以得到, 所对应的无向图是连通的. 因此, 假设3成立.
下面将使用三个实例来说明本文所提方法的有效性. 在所有的三个实例中, 算法1中具体的参数设置如下: 收集样本数据的实验次数$ \eta = 100>(3\;\times 4) / 2 = 6 $, 初始低增益参数$ \varepsilon = 1 $. 同时, 使用注5中的方法对$ \varepsilon $进行更新, 选择$ \alpha = 0.5 $. 后续将通过改变不同的初始状态来说明算法1的有效性.
例1. 在本例中, 将所有智能体的初始状态设置为$ [-1,\;1]\times [-1,\;1] $, 有限时域设置为$ \tau = 100 $, 然后将算法1应用于MASs (45)中, 最终得到低增益参数$ \varepsilon = 0.5 $. 同时, 将对应的最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $应用于系统中, 得到的6个智能体的系统状态如图2所示, 系统控制输入如图3所示. 从图2和图3可知, 通过算法1得到的控制输入可以实现有限时域一致性控制, 并且避免输入饱和.
例2. 在本例中, 将所有智能体的初始状态设置为$ [-10,\;10]\times [-10,\;10] $, 有限时域设置为$ \tau = 300 $, 然后将算法1应用于MASs (45)中, 最终得到低增益参数$ \varepsilon = 0.002 $. 同时, 将对应的最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $应用于系统中, 得到的6个智能体的系统状态如图4所示, 系统控制输入如图5所示. 不同于例1, 例2中智能体的初始状态的范围变大, 必然会影响MASs的一致性控制效果. 相比而言, 例2中智能体实现一致性控制的时间更长, 得到的低增益参数更小. 然而, 从图4和图5可知, 通过算法1得到的控制输入仍然可以实现有限时域一致性控制, 并避免输入饱和.
例3. 在本例中, 进一步加大了智能体初始状态的范围, 设置为$ [-100,\;100]\times [-100,\;100] $, 有限时域设置为$ \tau = 1\; 500 $, 然后将算法1应用于MASs (45)中, 最终得到低增益参数$ \varepsilon = 1.220\; 7\times 10^{-4} $. 同时, 将对应的最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $应用于系统中, 得到的6个智能体的系统状态如图6所示, 系统控制输入如图7所示. 从所得结果可知, 所提方法可以在有限时域内实现一致性控制, 并避免输入饱和.
以上三个例子证明了本文所提算法的有效性, 同时说明了如果智能体的初始状态越大, 控制输入需要配合越小的LGF控制增益矩阵$ K_{\varepsilon}(k) $以避免输入饱和, 因此低增益参数$ \varepsilon $将会迭代更多的次数, 从而得到更小的输入值. 此外, 在输入饱和度相等的情况下($ c = 1 $), 初始状态越大, 智能体实现一致性的速度越慢, 如图2、图4和图6所示. 通过以上三个例子, 也对定理3进行了验证.
3.2 仿真实验2
在本节将所提模型无关有限时域一致性控制算法与文献[38]针对执行器饱和的模型无关无限时域一致性控制方法进行对比.
考虑一个由5个智能体组成的离散时间MASs, 其动力学方程为(1), 相关的矩阵为:
$$ \begin{equation} A = \begin{bmatrix} 0.995 & -0.194\\ 0.194 & 0.995 \end{bmatrix},\;\quad B = \begin{bmatrix} 1 \\ 1 \end{bmatrix} \end{equation} $$ (46) 矩阵$ A $的特征值$ 0.980\; 1\pm0.198\; 7\mathrm{i} $都在单位圆上, 且$ (A,\; B) $是可控的. 因此, 假设2成立. 在本节仿真中, 执行器饱和函数的饱和阈值设为$ c = 1 $. 离散时间MASs的通信拓扑用图8所示的无向图表示. 从图8中可以得到, 所对应的无向图是连通的. 因此, 假设3成立.
图 11 例2中有限时域方法获得的一致性误差Fig. 11 Consensus errors obtained by finite-horizon method in example 2针对本文所提算法1的相关参数设置如下: 有限时域$ \tau = 120 $, 收集样本数据的实验次数$ \eta = 100 $, 初始低增益参数$ \varepsilon = 1 $, 低增益参数$ \varepsilon $调节参数$ \alpha = 0.9 $. 参考文献[38]所提无限时域算法的相关参数设置, 初始低增益参数$ \varepsilon = 1 $, $ M^0 = I $, $ K^0 = [0,\;0] $, 收集样本数据数量$ H = 100 $, 算法收敛参数设置为$ 0.000\; 01 $. 低增益参数$ \varepsilon $的更新规则和本文所提算法1一致.
例1. 在本例中, 首先设定5个智能体的初始状态为$ x_1(0) = [2.5,\, -2.5]^{\mathrm{T}} $, $ x_2(0) = [-1.5,\, 2]^{\mathrm{T}} $, $ x_3(0) = [-2,\; -3]^{\mathrm{T}} $, $ x_4(0) = [-2,\; -2]^{\mathrm{T}} $, $ x_5(0) = [1.5,\; 1.5]^{\mathrm{T}} $. 两种算法得到的最终低增益参数均为$ \varepsilon = 3.4\;\times 10^{-3} $. 采用文献[38]中所提算法得到的最优LGF控制增益矩阵$ K^*_{\varepsilon} = [-0.093\; 7,\; -0.073\; 0]^{\mathrm{T}} $. 将两种算法得到的最优LGF控制增益矩阵$ K_{\varepsilon}^*(k) $和$ K^*_{\varepsilon} $应用于MASs (23)中. 为了对比两种算法的一致性控制效果, 引入一致性控制误差$ \varepsilon_i(k) = \sum_{j = 1}^N d_{i j}(x_i(k)\;- x_j(k)) $. 仿真结果见图9和图10.
图 9 例1中有限时域方法获得的一致性误差Fig. 9 Consensus errors obtained by finite-horizon method in example 1
图 10 例1中无限时域方法获得的一致性误差Fig. 10 Consensus errors obtained by infinite-horizon method in example 1例2. 在本例中, 改变5个智能体的初始状态为$ x_1(0) = [1,\; 2]^{\mathrm{T}} $, $ x_2(0) = [-0.5,\; -0.1]^{\mathrm{T}} $, $ x_3(0) = [0.3, 2]^{\mathrm{T}} $, $ x_4(0) = [0.8,\; 0.2]^{\mathrm{T}} $, $ x_5(0) = [-3,\; -2]^{\mathrm{T}} $. 两种算法得到的最终低增益参数均为$ \varepsilon = 7.1\times 10^{-3} $. 文献[38]所提算法得到的最优LGF控制增益矩阵为$ K_{\varepsilon} = [-0.132\; 4,\; -0.110\; 6]^{\mathrm{T}} $. 最终得到的仿真结果见图11和图12.
图 12 例2中无限时域方法获得的一致性误差Fig. 12 Consensus errors obtained by infinite-horizon method in example 2另外, 本文用每个智能体对应一致性误差的绝对误差积分(Integral absolute error, IAE)的平均值和均方误差(Mean square error, MSE)的和两个指标[47−48]来评价本仿真实验的控制效果, 结果见表1.
表 1 对比实验评价指标Table 1 Evaluation indices of comparison experiments$100\le k \le 120$ ${\mathrm{IAE}}$ ${\mathrm{MSE}}$ 例1−有限时域方法 0.637 7 0.005 4 例1−无限时域方法 10.264 9 2.116 9 例2−有限时域方法 1.074 8 0.014 7 例2−无限时域方法 5.186 9 0.510 9 $$ \begin{equation} {\mathrm{IAE}} = \frac{\sum\limits_{i = 1}^N\sum\limits_{k = 0}^{k^*}|\varepsilon_i(k)|}{N} \nonumber \end{equation} $$ $$ \begin{equation} {\mathrm{MSE}} = \sum\limits_{i = 1}^N\frac{1}{k^*} \sum\limits_{k = 0}^{k^*}|\varepsilon_i(k)|^2\nonumber \end{equation} $$ 同时, 为了对比两种算法的一致性控制效果, 统计了智能体一致性误差对应的调节时间指标(以一致性误差范围的$ \pm\ 2\% $进行计算), 在不同初始状态下, 将时域参数均设置为200, 每个智能体对应的调节时间如表2和表3所示.
表 2 例1中一致性误差调节时间Table 2 Consensus error setting time in example 1例1−调节时间 有限时域方法 无限时域方法 智能体1 109 137 智能体2 119 161 智能体3 104 127 智能体4 109 137 智能体5 90 110 表 3 例2中一致性误差调节时间Table 3 Consensus error setting time in example 2例2−调节时间 有限时域方法 无限时域方法 智能体1 108 131 智能体2 116 158 智能体3 120 183 智能体4 108 131 智能体5 84 93 由图9 ~ 图12以及表1可知, 本文所提算法能够更快地实现一致性控制, 一致性误差较小. 同时由表2和表3可知, 在一定的时间范围内, 本文所提的有限时域一致性控制算法得到的一致性性能指标较文献[38]所提无限时域一致性控制算法要好, 这也说明了本文提出算法的优越性.
4. 结束语
本文提出一种基于Q学习的数据驱动算法, 用于求解具有未知模型参数、执行器饱和的离散时间MASs的有限时域一致性控制问题. 首先结合LGF方法, 将执行器饱和的有限时域一致性控制问题转化为执行器饱和的单智能体最优控制问题, 给出原问题的控制器设计方案. 然后在未知系统模型参数的条件下, 设计基于Q学习的数据驱动后向时间算法逼近求解MTVRE, 用以获取LGF控制增益矩阵, 并给出该算法的收敛性说明. 最后, 给出仿真结果来验证基于Q学习的有限时域一致性控制算法的有效性, 并证明智能体的初始状态会影响收敛速度的问题. 同时, 还给出对比实验来评价有限时域一致性控制算法与无限时域一致性控制算法的控制效果.
在本文提出的方法中, 有限时域参数 $ \tau $ 作为算法1的输入参数, 其在参数选择过程中需凭借经验来进行设定. 在未来的研究中, 将探讨更为精确的有限时域参数设置方法, 以确定 $ \tau $ 的边界条件, 从而设定合理的有限时域参数 $ \tau $.
-
图 3 3车场K-means未平衡聚类与平衡聚类比较
Fig. 3 Comparison of unbalanced K-means cluster and balanced K-means cluster of three depots
表 1 符号及定义
Table 1 Symbols and definitions
符号 释义 符号 释义 $F_1$ 运输距离费用 $H_{PM}$ 车场P 中有$H_{PM}$辆$M$类型的车辆 $F_2$ 车辆固定成本 $r(A)$ 完成客户子集$A$中所有客户的配送需要的最少车辆数 $F_3$ 燃油消耗费用 $N$ 总共有$N$个客户 $F_4$ 时间窗惩罚费用 $V$ 客户编号集合$\{1,2,\cdots,\ N\} $ (0 表示车场) $C_{M1}$ 第$M$种类型车辆的距离费用系数 $M_t$ 共有$M_t$种类型的车辆 $C_{M2}$ 第$M$种类型车辆的固定发车费用系数 $k$ 车辆编号 $C_{M3}$ 第$M$种类型车辆的燃油费用系数 $x_{PMijk}$ 车场$P$车型$M$的第$k$辆车从客户$i$到客户$j$的决策变量 $C_1$ 配送车辆提前到达的单位惩罚费用 $d_{ij}$ 客户$i$到客户$j$的距离 $C_2$ 配送车辆迟到的单位惩罚费用 ${ET}_i$ 客户$i$要求的最早到达时间 $i$ 客户点$i$ ${LT}_i$ 客户$i$要求的最晚到达时间 $j$ 客户点$j$ $S_i$ 客户$i$要求的卸货时间 $P$ $\{1,2,\cdots,\ P_t\} $ 车场编号$q_i$ 客户i 要求的货物需求量 $P_s$ 全部车场集合 $t_i$ 车辆到达客户$i$的时间 $P_t$ 总共有$P_t$个车场$P$ $M$ 车型编号 $M_s$ 全部车型集合$\{1,2,\cdots,\ M_t\}$ $Q_M$ 第$M$种车型的最大载重量 $H_{PMS}$ 车场$P$中车型$M$的全部车辆集合$\{1,2,\cdots,\ H_{PM}\}$ ${FU}_{Mij}$ 车型为$M$的车辆从客户$i$到客户$j$之间的耗油量 注: 综合燃油消耗模型中的其他相关参数设定参考文献 [25]. 
下载: 导出CSV
表 2 目标函数中的相关系数
Table 2 Coefficients in the object function
符号 数值 $C_{M1}$ 1.5 (元/km) $C_{M2}$ 300 ~ 800 (元/辆) $C_{M3}$ 7.6 (元/l) $C_{1}$ 15 (元/h) $C_{2}$ 20 (元/h)
下载: 导出CSV
表 3 主要参数与水平
Table 3 Main parameters and level
主要参数 水平设置 1 2 3 4 $\alpha$ 1.25 1.5 1.75 2.0 $\beta$ 10 1.5 2.0 2.5 $P_m$ 1.1 1.2 1.3 1.4 $W$ 500 1000 1500 2000
下载: 导出CSV
表 4 参数设置的正交表
Table 4 Orthogonal table of parameter settings
组合编号 水平设置 AVR (元) $\alpha$ $\beta$ $P_m$ $W$ 1 1 1 1 1 9677 2 1 2 2 2 9625 3 1 3 3 3 9613 4 1 4 4 4 9541 5 2 1 2 3 9745 6 2 2 1 4 9624 7 2 3 4 1 9602 8 2 4 3 2 9593 9 3 1 3 4 9836 10 3 2 4 3 9703 11 3 3 1 2 9654 12 3 4 2 1 9612 13 4 1 4 2 9865 14 4 2 3 1 9689 15 4 3 2 4 9656 16 4 4 1 3 9672
下载: 导出CSV
表 5 各参数不同水平下的平均响应值和影响力
Table 5 Average response values and influences table at different levels of each parameter
水平 水平设置 $\alpha$ $\beta$ $P_m$ $W$ 1 9614 9780 9656 9645 2 9641 9660 9659 9684 3 9701 9631 9683 9683 4 9720 9604 9677 9664 极差 106 176 27 39 影响力排名 2 1 4 3
下载: 导出CSV
表 6 4种不同车型相关参数设置
Table 6 Related parameter settings for four different vehicle types
车型
列表车型参数 载重量 (kg) 空车重量(kg) 平均速度(km/h) 固定费用(元) 最大承载货物数 (件) Type 1 200 1600 60 ~ 80 300 ~ 400 20 Type 2 500 2700 50 ~ 70 400 ~ 500 30 Type 3 600 3500 40 ~ 60 500 ~ 600 40 Type 4 800 5000 30 ~ 50 600 ~ 800 50
下载: 导出CSV
表 7 EACO_IBKA与其他算法对比结果
Table 7 Comparison results of EACO_IBKA with other algorithms
N_Pt EACO_IBKA EACO_KM EACO_NNA EACO1 DHACO ${ {T} }({{\rm{s}}} )$ 最优 平均 最差 标准差 最优 平均 最差 标准差 最优 平均 最差 标准差 最优 平均 最差 标准差 最优 平均 最差 标准差 48_2 11118 11800 12163 95 10745 11181 11579 97 11558 12080 12539 99 9650 10390 11026 87 10255 11068 11663 90 10 96_2 17483 18155 18549 183 18037 18379 19146 187 16768 17675 18229 190 15371 17009 17771 169 16011 17559 18161 174 19 144_2 24628 25435 26983 308 24880 25369 27201 314 25435 26710 27702 320 24366 25969 27419 318 24475 26884 28356 328 29 192_2 27522 28546 29649 411 28457 29482 30261 419 28758 29560 31019 428 28379 29838 30618 432 28524 30863 31665 445 38 240_2 31505 32699 34166 508 32517 33677 35235 518 32676 33723 35179 529 32722 34769 35906 534 33643 36475 37487 550 48 288_2 39217 41028 43326 592 40412 42363 44162 604 41179 42142 44696 616 41283 43164 44325 622 42606 43930 45137 641 58 360_2 53748 56268 58544 847 54672 57402 60199 864 55251 57847 59619 881 56276 58653 61938 890 56965 59409 62929 916 72 48_3 10767 11035 11572 83 10221 10827 11115 84 10995 11492 11929 93 9178 9883 10489 81 9754 10529 11095 83 14 96_3 16638 17278 17654 174 17066 17491 18222 182 15957 16821 17349 180 14627 16582 17912 191 15236 17186 18141 184 29 144_3 23443 24211 25685 293 23659 24149 25893 290 24211 25426 26371 302 23194 24720 26101 296 23297 25592 26993 305 43 192_3 26199 27175 28225 392 27090 28066 28826 399 27376 28140 29529 403 27016 28405 29147 407 27153 29381 30145 420 58 240_3 29992 31130 32527 484 30956 32061 33545 494 31108 32105 33491 499 31151 33101 34184 504 32029 34726 35690 519 72 288_3 37337 39062 41251 564 38476 40333 42047 575 39205 41123 42555 581 39305 41096 42202 587 40565 41826 42975 604 86 360_3 51176 53576 55744 806 52056 54656 57320 823 52608 55080 56768 831 53584 55848 58976 839 54240 56568 59920 864 108 48_4 10239 10617 11313 82 9943 10418 10886 84 10807 11144 11671 86 8939 9615 9883 78 9545 10051 10539 80 19 96_4 16062 16943 17257 166 16810 17146 17481 175 15871 16425 17434 160 16012 16907 17336 179 16391 16991 17496 177 38 144_4 22700 23723 25488 279 23650 24242 25187 291 23765 24990 25945 296 22747 24273 25602 288 22841 25114 26506 294 58 192_4 25690 26612 27737 373 26528 27235 27938 383 26846 27620 28988 391 26443 27896 28628 403 26623 28850 29614 411 77 240_4 29375 30447 31877 461 30317 31401 32874 473 30447 31411 32863 483 30545 32484 33502 511 31411 34401 34997 521 96 288_4 36519 38266 40433 537 37624 39493 41240 549 38365 39338 41737 560 38442 40278 41339 583 39692 40974 42135 595 115 360_4 50480 52904 55056 768 51344 53960 56608 786 51904 54392 56088 802 52088 55184 58312 833 53584 55904 59216 850 144 平均值 28183 29377 30724 400 28831 29968 31284 409 29100 30250 31510 416 28634 30289 31553 421 29278 31156 32422 431 —
下载: 导出CSV
表 8 HKMA与其他划分算法的对比结果
Table 8 Comparison results of HKMA and the other dividing algorithms
$ {{N}}\_ {{M}}_{{t}}$ EACO_HKMA EACO_RDA EACO_RAA EACO_KEW EACO2 TSA_RDA $ {T}({{{\rm{s}}}})$ 最优 平均 最差 最优 平均 最差 最优 平均 最差 最优 平均 最差 最优 平均 最差 最优 平均 最差 48_2 12801 13612 13901 13220 14646 15128 12467 13293 13998 12607 13317 13922 12218 13798 14202 13352 14792 15279 15 96_2 16299 16759 17342 17570 19461 20106 16238 17661 18591 16754 17704 18101 16543 17838 18777 17746 19656 20307 29 144_2 22061 23265 24089 23361 25860 26717 22390 23665 24721 22266 23533 24590 22838 23934 25215 23595 26119 26984 44 192_2 24847 25998 26933 25983 26905 28228 26000 27201 28372 26039 26918 28238 26520 27745 28939 26243 27174 28510 57 240_2 25958 27124 28797 26951 28002 29726 29253 30130 30753 29409 30076 30878 30131 31034 31976 27221 28282 30023 72 288_2 32225 33356 33838 33783 34713 36290 32740 35242 37189 32540 35156 37175 34050 36652 38677 34121 35060 36653 87 360_2 48344 50050 51763 50695 52088 54440 49129 52874 55793 48828 52743 53780 51094 54989 58025 51202 52609 54984 108 48_3 11896 12400 13111 12574 13919 14386 11755 11999 13312 11995 12658 13236 11520 12240 13999 12700 14058 14530 21 96_3 15777 15929 16997 16704 18494 19122 16011 16792 17680 15932 16839 17594 16171 16960 17857 16871 18679 19313 44 144_3 20965 21171 22661 22207 24587 25394 21291 22312 23490 21170 22366 23380 21717 22758 23960 22429 24833 25648 65 192_3 23624 24718 25597 24704 25570 26823 24714 25850 26974 24748 25588 26836 25455 26626 27783 24951 25826 27091 87 240_3 24680 25776 27367 25620 26616 28248 27804 28634 29235 27948 28582 29351 28916 29779 30404 25876 26882 28530 108 288_3 30621 31700 32158 32106 32984 34483 31111 33493 35335 30929 33414 35335 32667 35168 37102 32427 33314 34828 129 360_3 45933 47554 48245 48172 49502 51726 46691 50237 53019 47401 50118 51004 49026 52749 55670 48654 49997 52243 162 48_4 10755 10998 12257 11330 12533 12956 10700 11393 11985 10800 11402 11927 10999 11766 11999 11443 12658 13086 29 96_4 14214 14349 14672 15047 16655 17221 14422 15118 15927 14346 15171 15843 14566 15269 16086 15197 16822 17393 57 144_4 18878 19057 20506 19998 22039 22865 19165 20085 21148 19057 20140 21048 19548 20487 21571 20198 22360 23094 87 192_4 21273 22251 23043 22237 23018 24145 22257 23268 24281 22282 23034 24158 22925 23966 25009 22459 23248 24386 116 240_4 22215 23209 24642 23073 23966 25438 25035 25778 26316 25161 25735 26431 26036 26809 27369 23304 24206 25692 144 288_4 27563 28542 29950 28902 29689 31040 28012 30149 31812 27842 30088 31815 29413 31656 33403 29191 29986 31350 173 360_4 41350 42809 44034 43360 44558 46561 42028 45218 47722 41776 43111 45711 44129 47479 50108 43794 45004 47027 216 平均值 24407 25197 25970 25600 26948 28145 25201 26659 27983 25230 26652 27636 26023 27605 28944 25856 27217 28426 —
下载: 导出CSV
表 9 EACO_CD性能验证
Table 9 Performance verification of EACO_CD
$ {{N}}\_{{ P}}_{{t}}\_{{M}}_{{t}}$ EACO_CD (EACO_IBKA_HKMA) IACO_CD (IACO_NNA_SWA) IHGA ${ {T} }({ {{\rm{s}}} })$ 最优值 平均值 最差值 标准差 最优值 平均值 最差值 标准差 最优值 平均值 最差值 标准差 48_2_2 12145 12478 12859 232 12880 13212 13509 241 11632 11967 12227 220 15 96_2_2 16370 16526 16904 295 17800 18071 18381 272 15852 16011 16376 285 31 144_2_2 21010 21305 21704 324 24589 24825 25406 336 21049 21640 22154 348 51 192_2_2 23664 24760 25650 556 26040 27250 28218 576 24863 26010 26948 590 71 240_2_2 24722 25832 27426 588 26707 27903 29640 641 27444 28720 29947 692 96 288_2_2 30690 31768 32227 724 33772 34962 35457 803 35311 36550 37076 824 123 360_2_2 39903 41303 41906 951 43907 45448 46111 1059 45507 46098 47781 1197 162 48_3_2 12999 14064 14425 293 13656 14787 15372 336 12742 13802 14337 304 21 96_3_2 16883 17869 18993 358 18568 20115 20999 343 16800 18283 19011 363 48 144_3_2 21916 23728 24668 456 25647 27774 28874 539 23029 24921 25906 491 76 192_3_2 24257 25492 26408 668 26692 28059 29065 671 26688 28052 29053 746 109 240_3_2 24672 26558 28026 792 27884 30021 31675 740 28382 30554 32236 787 144 288_3_2 30258 31572 32860 907 34814 36319 37797 1017 35414 36956 38450 1092 183 360_3_2 39347 41061 42734 1141 43291 45187 47020 1263 47228 49291 51297 1389 243 48_3_3 12284 13225 13755 280 12830 13906 14446 322 12123 13253 13797 287 29 96_3_3 15867 17082 17860 388 17470 18911 19664 368 16514 17874 18584 478 62 144_3_3 20620 22321 23193 438 22488 24336 25294 490 24128 26136 27154 522 102 192_3_3 22816 23977 24832 631 25108 26383 27320 708 26242 27588 28567 743 145 240_3_3 24209 25982 27348 655 26693 28731 30318 770 27859 29995 31632 803 192 288_3_3 28458 29693 30892 828 32742 34157 35530 957 35582 37120 38622 1048 245 360_3_3 37006 38610 40166 1079 40722 42474 44186 1198 47081 49427 51523 1399 288 平均值 23814 25010 25945 599 26395 27754 28775 650 26737 28107 29175 696 —
下载: 导出CSV
-
[1] Dantzig G B, Ramser J H. The truck dispatching problem. Management Science, 1959, 6(1): 80−91 doi: 10.1287/mnsc.6.1.80 [2] Anbuudayasankar S P, Ganesh K, Mohapatra S. Models for Practical Routing Problems in Logistics. Cham: Springer, 2014. [3] Li H Q, Yuan J L, Lv T, Chang X Y. The two-echelon time-constrained vehicle routing problem in linehaul-delivery systems considering carbon dioxide emissions. Transportation Research Part D: Transport and Environment, 2016, 49: 231−245 doi: 10.1016/j.trd.2016.10.002 [4] 赵燕伟, 张景玲, 王万良. 物流配送的车辆路径优化方法. 北京: 科学出版社, 2014.Zhao Yan-Wei, Zhang Jing-Ling, Wang Wan-Liang. Vehicle Routing Optimization Methods for Logistics Distribution. Beijing: Science Press, 2014. [5] Sbihi A, Eglese R W. Combinatorial optimization and green logistics. 4OR, 2007, 5(2): 99−116 doi: 10.1007/s10288-007-0047-3 [6] Chen D S, Batson R G, Dang Y. Applied Integer Programming. Hoboken: John Wiley & Sons, 2010. [7] Jabir E, Panicker V V, Sridharan R. Design and development of a hybrid ant colony-variable neighbourhood search algorithm for a multi-depot green vehicle routing problem. Transportation Research Part D: Transport and Environment, 2017, 57: 422−457 doi: 10.1016/j.trd.2017.09.003 [8] Kaabachi I, Jriji D, Krichen S. An improved ant colony optimization for green multi-depot vehicle routing problem with time windows. In: Proceedings of the 18th IEEE/ACIS International Conference on Software Engineering, Artificial Intelligence, Networking and Parallel/Distributed Computing (SNPD). Kanazawa, Japan: IEEE, 2017. [9] Xiao Y Y, Konak A. The heterogeneous green vehicle routing and scheduling problem with time-varying traffic congestion. Transportation Research Part E: Logistics and Transportation Review, 2016, 88: 146−166 doi: 10.1016/j.tre.2016.01.011 [10] Kwon Y J, Choi Y J, Lee D H. Heterogeneous fixed fleet vehicle routing considering carbon emission. Transportation Research Part D: Transport and Environment, 2013, 23: 81−89 doi: 10.1016/j.trd.2013.04.001 [11] Geetha S, Vanathi P T, Poonthalir G. Metaheuristic approach for the multi-depot vehicle routing problem. Applied Artificial Intelligence, 2012, 26(9): 878−901 doi: 10.1080/08839514.2012.727344 [12] Geetha S, Poonthalir G, Vanathi P T. Nested particle swarm optimisation for multi-depot vehicle routing problem. International Journal of Operational Research, 2013, 16(3): 329−348 doi: 10.1504/IJOR.2013.052336 [13] Ho W, Ho G T S, Ji P, Lau H C W. A hybrid genetic algorithm for the multi-depot vehicle routing problem. Engineering Applications of Artificial Intelligence, 2008, 21(4): 548−557 doi: 10.1016/j.engappai.2007.06.001 [14] Wang Y, Assogba K, Liu Y, Ma X L, Xu M Z, Wang Y H. Two-echelon location-routing optimization with time windows based on customer clustering. Expert Systems With Applications, 2018, 104: 244−260 doi: 10.1016/j.eswa.2018.03.018 [15] Dondo R, Cerdá J. A cluster-based optimization approach for the multi-depot heterogeneous fleet vehicle routing problem with time windows. European Journal of Operational Research, 2007, 176(3): 1478−1507 doi: 10.1016/j.ejor.2004.07.077 [16] Tang Y L, Cai Y G, Yang Q J. Improved ant colony optimization for multi-depot heterogeneous vehicle routing problem with soft time windows. Journal of Southeast University (English Edition), 2015, 31(1): 94−99 [17] Dorigo M, Maniezzo V, Colorni A. Ant system: Optimization by a colony of cooperating agents. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics), 1996, 26(1): 29−41 doi: 10.1109/3477.484436 [18] 王素欣, 高利, 崔小光, 曹宏美. 多需求点车辆调度模型及其群体智能混合求解. 自动化学报, 2008, 34(1): 102−104Wang Su-Xin, Gao Li, Cui Xiao-Guang, Cao Hong-Mei. Study on multi-requirement points vehicle scheduling model and its swarm mix algorithm. Acta Automatica Sinica, 2008, 34(1): 102−104 [19] Lee C Y, Lee Z J, Lin S W, Ying K C. An enhanced ant colony optimization (EACO) applied to capacitated vehicle routing problem. Applied Intelligence, 2010, 32(1): 88−95 doi: 10.1007/s10489-008-0136-9 [20] Yu B, Yang Z Z. An ant colony optimization model: The period vehicle routing problem with time windows. Transportation Research Part E: Logistics and Transportation Review, 2011, 47(2): 166−181 doi: 10.1016/j.tre.2010.09.010 [21] Ding Q L, Hu X P, Sun L J, Wang Y Z. An improved ant colony optimization and its application to vehicle routing problem with time windows. Neurocomputing, 2012, 98: 101−107 doi: 10.1016/j.neucom.2011.09.040 [22] Yan F L. Autonomous vehicle routing problem solution based on artificial potential field with parallel ant colony optimization (ACO) algorithm. Pattern Recognition Letters, 2018, 116: 195−199 doi: 10.1016/j.patrec.2018.10.015 [23] 陈希琼, 胡大伟, 杨倩倩, 胡卉, 高扬. 多目标同时取送货车辆路径问题的改进蚁群算法. 控制理论与应用, 2018, 35(9): 1347−1356 doi: 10.7641/CTA.2018.80085Chen Xi-Qiong, Hu Da-Wei, Yang Qian-Qian, Hu Hui, Gao Yang. An improved ant colony algorithm for multi-objective vehicle routing problem with simultaneous pickup and delivery. Control Theory & Applications, 2018, 35(9): 1347−1356 doi: 10.7641/CTA.2018.80085 [24] Xu H T, Pu P, Duan F. Dynamic vehicle routing problems with enhanced ant colony optimization. Discrete Dynamics in Nature and Society, 2018, 2018: Article No. 1295485 [25] Demir E, Bektaš T, Laporte G. An adaptive large neighborhood search heuristic for the pollution-routing problem. European Journal of Operational Research, 2012, 223(2): 346−359 doi: 10.1016/j.ejor.2012.06.044 [26] Bektaš T, Gouveia L. Requiem for the Miller-Tucker-Zemlin subtour elimination constraints. European Journal of Operational Research, 2014, 236(3): 820−832 doi: 10.1016/j.ejor.2013.07.038 [27] Toth P, Vigo D. 车辆路径问题. 北京: 清华大学出版社, 2011.Toth P, Vigo D. The Vehicle Routing Problem. Beijing: Tsinghua University Press, 2011. [28] Wolsey L A. Integer Programming. New York: Wiley, 1998. [29] Schneider M, Stenger A, Goeke D. The electric vehicle-routing problem with time windows and recharging stations. Transportation Science, 2014, 48(4): 500−520 doi: 10.1287/trsc.2013.0490 [30] Gu S S. A polynomial time solvable algorithm to linearly constrained binary quadratic programming problems with Q being a tri-diagonal matrix. In: Proceedings of the 5th International Conference on Intelligent Control and Information Processing (ICCIP). Dalian, China: IEEE, 2011. [31] Meng X H, Li J, Zhou M C, Dai X Z, Dou J P. Population-based incremental learning algorithm for a serial colored traveling salesman problem. IEEE Transactions on Systems, Man, and Cybernetics: Systems, 2018, 48(2): 277−288 doi: 10.1109/TSMC.2016.2591267 [32] Wang Y, Ma X L, Lao Y T, Wang Y H. A fuzzy-based customer clustering approach with hierarchical structure for logistics network optimization. Expert Systems With Applications, 2014, 41(2): 521−534 doi: 10.1016/j.eswa.2013.07.078 [33] Wang Y, Zhang J, Assogba K, Liu Y, Xu M Z, Wang Y H. Collaboration and transportation resource sharing in multiple centers vehicle routing optimization with delivery and pickup. Knowledge-Based Systems, 2018, 160: 296−310 doi: 10.1016/j.knosys.2018.07.024 [34] Wang B D, Miao Y W, Zhao H Y, Jin J, Chen Y Z. A biclustering-based method for market segmentation using customer pain points. Engineering Applications of Artificial Intelligence, 2016, 47: 101−109 doi: 10.1016/j.engappai.2015.06.005 [35] Ji J C, Pang W, Zhou C G, Han X, Wang Z. A fuzzy k-prototype clustering algorithm for mixed numeric and categorical data. Knowledge-Based Systems, 2012, 30: 129−135 doi: 10.1016/j.knosys.2012.01.006 [36] Wang Y, Ma X L, Liu M W, Gong K, Liu Y, Xu M Z, et al. Cooperation and profit allocation in two-echelon logistics joint distribution network optimization. Applied Soft Computing, 2017, 56: 143−157 doi: 10.1016/j.asoc.2017.02.025 [37] He R H, Xu W B, Sun J X, Zu B Q. Balanced k-means algorithm for partitioning areas in large-scale vehicle routing problem. In: Proceedings of the 3rd International Symposium on Intelligent Information Technology Application. Nanchang, China: IEEE, 2009. [38] 王东风, 孟丽. 粒子群优化算法的性能分析和参数选择. 自动化学报, 2016, 42(10): 1552−1561 doi: 10.16383/j.aas.2016.c150774Wang Dong-Feng, Meng Li. Performance analysis and parameter selection of PSO algorithms. Acta Automatica Sinica, 2016, 42(10): 1552−1561 doi: 10.16383/j.aas.2016.c150774 [39] Beasley J E. Route first―Cluster second methods for vehicle routing. Omega, 1983, 11(4): 403−408 doi: 10.1016/0305-0483(83)90033-6 [40] Gillett B E, Miller L R. A heuristic algorithm for the vehicle-dispatch problem. Operations Research, 1974, 22(2): 340−349 doi: 10.1287/opre.22.2.340 [41] Yu B, Yang Z Z, Yao B Z. An improved ant colony optimization for vehicle routing problem. European Journal of Operational Research, 2009, 196(1): 171−176 doi: 10.1016/j.ejor.2008.02.028 [42] Mladenović N, Hansen P. Variable neighborhood search. Computers & Operations Research, 1997, 24(11): 1097−1100 [43] 李进, 傅培华. 具有固定车辆数的多车型低碳路径问题及算法. 计算机集成制造系统, 2013, 19(6): 1351−1362 doi: 10.13196/j.cims.2013.06.189.lij.007Li Jin, Fu Pei-Hua. Heterogeneous fixed fleet low-carbon routing problem and algorithm. Computer Integrated Manufacturing Systems, 2013, 19(6): 1351−1362 doi: 10.13196/j.cims.2013.06.189.lij.007 期刊类型引用(2)
1. 谭福容，孙绍伦，张森，陈先中，赵宝永. 基于泊松算法和多尺度特征编码网络的三维料面重构及修复. 冶金自动化. 2024(02): 94-102 . 
百度学术2. 赵炯. 高炉炉顶气密箱布料溜槽不同工况下的使用特性分析及优化. 山西冶金. 2024(09): 124-125+128 . 百度学术其他类型引用(2)
-
下载:
计量
- 文章访问数: 2510
- HTML全文浏览量: 281
- PDF下载量: 249
- 被引次数: 4
