[插电式混合动力车辆][交替方向乘子法（ADMM）结合CVX]插电式混合动力车辆的能源管理：基于凸优化算法用于模型预测控制MPC研究（Matlab代码实现）

 💥💥💞💞欢迎来到本博客❤️❤️💥💥

🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者。

⛳️座右铭：行百里者，半于九十。

📋📋📋本文目录如下：🎁🎁🎁

目录

 ⛳️赠与读者

💥1 概述

📚2 运行结果

🎉3 参考文献

🌈4 Matlab代码、数据、文章

---

 ⛳️赠与读者

👨‍💻做科研，涉及到一个深在的思想系统，需要科研者逻辑缜密，踏实认真，但是不能只是努力，很多时候借力比努力更重要，然后还要有仰望星空的创新点和启发点。当哲学课上老师问你什么是科学，什么是电的时候，不要觉得这些问题搞笑。哲学是科学之母，哲学就是追究终极问题，寻找那些不言自明只有小孩子会问的但是你却回答不出来的问题。建议读者按目录次序逐一浏览，免得骤然跌入幽暗的迷宫找不到来时的路，它不足为你揭示全部问题的答案，但若能让人胸中升起一朵朵疑云，也未尝不会酿成晚霞斑斓的别一番景致，万一它居然给你带来了一场精神世界的苦雨，那就借机洗刷一下原来存放在那儿的“躺平”上的尘埃吧。

     或许，雨过云收，神驰的天地更清朗.......🔎🔎🔎

💥1 概述

摘要：
本文详细介绍了对用于求解凸优化问题的算法在混合动力车辆能源管理中的计算性能进行研究，其中包括非线性损耗。提出了一种投影内点法，通过在控制输入上应用不等式约束作为投影，减少了牛顿步骤矩阵求逆的大小和复杂性，并通过与交替方向乘子法（ADMM）算法和通用凸优化软件CVX的比较在仿真中展示了其性质。发现当需要适度精度的解时，ADMM算法具有有利的性质，而当需要高精度时，投影内点法更有利，并且两者都比CVX快得多。

增加电动化的道路车辆被确定为解决重要社会问题的关键短期方案，包括气候变化和空气污染[1]。插电式混合动力电动车（PHEV），其中电动推进系统与内燃发动机相辅相成，目前是一种常见配置。尽管锂离子电池的低能量密度和长充电时间限制了全电动动力总成的可行性，但对日常驾驶行为的分析表明，使用纯电动续航里程仅为40英里的混合动力车辆可以用电驱动50%的内燃机里程[2]。然而，额外动力源的引入带来了一个棘手的问题：在给定旅程的每一刻，应该从电动机提供多少功率，从发动机提供多少功率？

这被称为能量管理问题[3]，一个简单的启发式方法是充电耗尽/充电维持策略，即只有在电池充分耗尽时才从电动机提供动力，然后在旅程结束前以充电维持模式运行车辆[2]。然而，已经证明通过同时从电动机和发动机提供动力，并在整个旅程中调节从每个动力源提供的比例，可以实现燃油消耗的显著节省，这被称为“混合模式”[4]。有几种控制动力总成的方法，可以通过动态规划获得全局最优解，用于复杂的非线性模型和/或整数控制决策（如换挡选择）[5]，[6]，但这种方法对于在线解决方案来说计算要求过高。其他研究已经调查了基于庞特里亚金最小原理的方法[7]–​[9]，尽管在仍保证最优性的情况下，强制执行复杂约束（如总电荷状态约束或发动机切换）是具有挑战性的。

模型预测控制（MPC）在这个应用中显示出潜力，因为它对车辆模型和未来驾驶行为预测的不确定性具有固有的鲁棒性[10]，并且整个混合动力总成中可以使用非线性损耗模型以提高性能[11]。当考虑到换挡和发动机切换时，相关的在线优化问题仍然难以计算，因此这些元素通常会被从问题中移除或在外部进行优化[11]–​[14]，并且已经证明仅通过制定功率平衡，使用非线性损失而无需简化，可以将问题表述为具有线性状态动力学的凸优化问题[15]。

在过去的30年里，解决不等式约束优化问题最流行的算法一直是内点法[16]。最初被构造为“原始”方法，通过用对数障碍函数逼近优化问题中的不等式约束，它们固有的病态和数值效率低下使内点法在卡马卡尔（Karmakar）于1984年发表其方法[17]之前无效。随后发展出的“原始-对偶”内点法[18]表现出优秀的理论和实际特性，包括多项式复杂度和在问题规模变化时迭代次数近似恒定，如今，在MPC优化领域，大量研究成果致力于原始-对偶内点法[19]–[21]的发展和应用。

在最近发表的文献中，尚未研究解决能源管理问题的凸优化形式的算法，通常使用通用凸优化软件[13]，[14]，[22]，[23]来解决优化问题，除了在[15]中提出的交替方向乘子法（ADMM）算法。本文的第一个贡献是针对这个问题提出了一个投影内点求解器，通过将决策变量上的逐元素不等式约束作为投影，降低了与原始-对偶牛顿步骤相关的矩阵求逆的规模，从而降低了每次迭代的计算需求。该算法不特定于领域，适用于具有线性动力学、控制变量的可分离凸成本函数以及控制和状态变量的上下界的任何MPC优化问题。

本文的一个主要动机是确定对于凸PHEV能源管理制定，二阶和一阶方法的相对计算优势，第二个贡献是一组数值研究，比较了投影内点算法与[15]中的ADMM算法的性能。在这些研究中，我们证明了投影内点方法具有优越的收敛性质，但需要更多时间才能获得具有适度精度的解决方案，因此，仅适用于较短时间范围内的实时解决方案（在本例中，少于500个样本）。我们还证明了两种方法都比CVX[24]快得多，而在ADMM中（使用[15]中的改进实现），我们展示了第一种能够在实时情况下解决能源管理问题的方法，在考虑非线性系统动力学并在整个时间范围内强制执行功率和电荷状态的硬限制时，可以在长时间范围（≥1000个样本）内实现。

本文的结构如下。第二节定义了能源管理问题、MPC框架和凸重构，第三节详细介绍了投影内点方法。第四节阐明了[15]中的ADMM算法，第五节展示了数值实验，第六节对本文进行了总结。

📚2 运行结果

通过运行Main.m脚本生成的图表如下所示，其中ADMM的结果与CVX的结果完全匹配。，ADMM代码大约需要0.6秒，而CVX则需要大约800秒。

部分代码：

function [ E, u, time, iterations ] = f_ADMM( coeffs,Pdrv,E0,Pbmin,Pbmax,xmin,xmax,P,C,R,V,misc )

N = length(Pbmin);

rho1 = 2.34E-4;
rho2 = 1E-8;

alpha2 = coeffs(:,1);
alpha1 = coeffs(:,2);
alpha0 = coeffs(:,3);
beta2 = coeffs(:,4);
beta1 = coeffs(:,5);
beta0 = coeffs(:,6);

%% Feasibility check

%% Initialization (these operations can be performed offline and are not timed)

u = zeros(N,1);
u(C) = Pbmin(C);
zeta = zeros(N,1);
lambda1 = zeros(N,1);
lambda2 = zeros(N,1);

I = eye(N);
Psi = tril(ones(N,N));
M = (rho1/rho2 * inv(Psi) * (inv(Psi)') + I);
L = chol(M)';
L = sparse(L);
Diff = sparse(inv(Psi));
Difft = sparse(inv(Psi)');

%% Algorithm
tic

iterations = 0;
flag = 1;

while flag
    
    % The u update for k in P is the solution of a convex optimization
    % problem. This is solved using a newton method that is included in the
    % function f_BacktrackingNewtonVector
    u(P) = f_BacktrackingNewtonVector(alpha0(P), alpha1(P), alpha2(P), beta0(P), beta1(P), beta2(P), V, R, Pdrv(P), rho1, zeta(P), lambda1(P), Pbmin(P), Pbmax(P));
    
    % The x update is trivial
    x = E0 - cumsum(zeta) - lambda2;
    x(x > xmax) = xmax;
    x(x < xmin) = xmin;
    
    %hold the current value of zeta for residual calculations
    zetahold = zeta;
    
    % The zeta update is solved using a method detailed in Appendix D of 
    % "Optimal Power Allocation in Battery/Supercapacitor Electric Vehicles 
    % using Convex Optimization", available at
    % https://ieeexplore.ieee.org/document/9193947
    % and
    % https://arxiv.org/abs/2005.03678
    vec = rho1 * (u + lambda1) - rho2 * cumsum(x - E0 + lambda2, 'reverse');
    vec = vec / rho2;
    vec = Difft * vec; 
    vec = Diff * vec; 
    vec = L \ vec;
    zeta = L' \ vec;
    
    % The residual and Lagrange multiplier updates are also trivial
    r = [u - zeta; x + cumsum(zeta) - E0];
    s = [rho1 * (zetahold - zeta); rho2*cumsum(zetahold - zeta)];
    lambda1 = lambda1 + (u - zeta);
    lambda2 = lambda2 + (x + cumsum(zeta) - E0);
    
    iterations = iterations + 1;
    
    % Termination criteria
    if iterations > misc.maxIterations
        flag = 0;
    end
    
    if max(norm(r), norm(s)) < misc.epsilon
        flag = 0;
    end
    
end

time = toc;

E = E0 - cumsum(u);

return
 

🎉3 参考文献

文章中一些内容引自网络，会注明出处或引用为参考文献，难免有未尽之处，如有不妥，请随时联系删除。

🌈4 Matlab代码、数据、文章