ARTICLE · INTELLIGENCE

战地情报 · 详情页

来自尧图项目组的一线实战观察与深度解析

基于NSGA-II的固定翼无人机混合电力推进系统多目标优化

基于NSGA-II的固定翼无人机混合电力推进系统多目标优化

搞过固定翼无人机混合电力推进系统的朋友,应该都体会过这种纠结:电池和燃油怎么配比,才能在满足航时航程的前提下,既让飞机尽量轻,又让油耗和电池损耗尽量低?这几个目标互相打架,单目标优化根本没法给答案。我最近基于NSGA-II算法做了一个模拟项目X,对一套混动固定翼推进系统做多目标优化,把完整Matlab代码跑通并整理了出来。这篇文章就把整个项目从建模、算法设计到代码实现、问题排查的思路全部拆开讲,给准备做类似飞行器能源系统优化的同学一个可以直接参考的框架。

1. 为什么固定翼混动系统优化必须用多目标算法

很多文章讲NSGA-II时,习惯先用ZDT或DTLZ这类测试函数演示一下帕累托前沿,贴一张图画风很漂亮,但离工程落地差得很远。工程问题的难点从来不在算法本身,而在于目标函数和系统模型的构建能否反映真实物理规律。飞行器混合电力推进系统恰好多目标特性极其明显。

我这个模拟项目X的推进架构是电池组和发电机并联,共同给直流母线供电,直流母线再经电机控制器驱动推进电机带动螺旋桨。系统里有几个核心参数需要决策:电池容量、发电机额定功率、电机额定功率。这三个参数直接决定三个互相矛盾的目标:

  • 起飞总质量要最小化,飞机越轻,结构成本和飞行性能越好。
  • 任务总油耗要最小化,燃油经济性决定航时和运营成本。
  • 电池寿命衰减系数要最小化,电池循环寿命直接影响全寿命周期成本。

如果给这三个目标加权重变成一个单目标函数,会遇到一个非常头疼的问题:你根本不知道权重设多少才算合理,而且目标量纲差异巨大,结果对权重极其敏感,换个权重答案就面目全非。更关键的是,单目标只能给你一个解,你得反复试不同权重才能知道设计空间中风险和性能的妥协关系。NSGA-II的意义在于它一次性输出一整条帕累托前沿,前沿上每个点都是一种不被其他配置完全压制的设计,决策者可以根据实际任务需求在这个解集里挑最合适的点。

我选用Matlab R2021b作为实现环境。虽然优化工具箱自带gamultiobj函数,本质也是NSGA-II的一种实现,我最终还是选择手写主流程。原因是工程优化里经常需要调整交叉变异算子、加入约束处理策略、甚至把决策变量映射到离散档位,内置函数在这些场景下并不顺手。手写代码的另一大好处是每一步都在掌控之中,结果不理想时能迅速定位是算子的问题还是目标函数的问题。

2. 混合电力推进系统的建模逻辑与关键能量耦合

做这类优化,最怕的就是把目标函数写成一个凭空捏造的漂亮公式,然后丢给算法跑,跑完发现结果毫无物理意义。混合电力推进系统的建模才是整个项目真正花时间的地方。

2.1 各部件模型与飞行任务剖面

系统模型我拆成几大块:飞行动力学、功率需求、电池模型、发电机模型。飞行动力学采用水平飞行模型,升力等于重力,推力等于阻力,这是经典计算方式。空气动力学阻力包含零升阻力和诱导阻力两部分,表达式为:

D = 0.5 * rho * V^2 * S * (CD0 + K * (L / (0.5 * rho * V^2 * S))^2)

电机输出功率由推力乘以速度再除以螺旋桨效率和电机效率得到。任务剖面我简化成三段:爬升段、巡航段、降落段。爬升段功率需求大,巡航段相对平稳。每个阶段的飞行参数和持续时间都不同,这些参数会直接影响系统对电池和发电机的需求分配。

能量管理策略是系统模型中的关键模块。我设定爬升阶段电池承担功率比例为0.8,巡航阶段比例降为0.3。这个设定有物理依据:爬升需要快速释放大功率,电池动态响应快且能短时高倍率放电;巡航阶段功率需求平缓,让发电机承担基础负荷更经济,电池只做削峰填谷。电池模型里我用了简化的循环老化模型,衰减系数与放电深度DOD和放电倍率C-rate挂钩,大电流放电会加速老化,这和现实中锂电池的规律是一致的。发电机模型相对简单,核心参数是燃油消耗率SFC和发电机效率,输出轴功率乘以比油耗即得燃油消耗。

2.2 质量-能量耦合迭代:最容易被忽略的坑

这个模型里有一个非常隐蔽但极重要的耦合关系:巡航段重量W_cruise里包含燃油重量m_fuel,而m_fuel由巡航段发电功率P_gen_cruise决定,P_gen_cruise由巡航功率P_cruise决定,P_cruise又依赖于W_cruise。这是个典型的闭环耦合,简化计算时如果只用第一轮估算的重量去算功率需求,不再迭代,总质量误差会非常可观。

我第一次跑通代码时没有做这个迭代,结果帕累托前沿上有几个点明显不合理,某些配置计算出的起飞质量看似很低,但仔细检查发现燃油重量没有反馈回巡航阻力计算,导致功率需求偏低、油耗估计偏小,整个质量估算链条不闭合。后来我在目标函数里加了固定点迭代:用当前总重量算功率需求,由功率需求算燃油消耗,再把燃油消耗加回总重量,重复2到3次。经过两次迭代后结果趋向稳定,前沿形状也合理了。

这个迭代细节建议所有做飞行器能量系统优化的同学重点关注。飞行器对重量极其敏感,忽略这个耦合循环会导致质量估算偏差10%以上,直接误导设计决策。地面车辆或船舶系统对重量没那么敏感,但这个思路依然适用。参数上我用的是一套2米翼展级别验证机的逆向参数,巡航速度30m/s、机翼面积1.2平方米、巡航时间3600秒。不同机型必须针对性地重新标定参数。

3. NSGA-II核心流程与Matlab实现细节

NSGA-II的核心逻辑是一个“进化加筛选”的循环:初始化种群,计算目标函数,非支配排序分层,计算拥挤度距离,锦标赛选择,交叉变异生成子代,合并父子种群后精英保留。算法本身不复杂,但每个环节的代码实现都有关键细节。

3.1 主循环与参数设定

种群规模我设为100,迭代200代,决策变量3个分别是电池容量、发电机额定功率、电机额定功率。随机种子固定为rng(2025),保证实验结果可复现。这个细节看上去简单,做实验管理时却非常重要,特别是调参过程中需要对比不同参数的影响,如果没有固定种子,每次结果的差异可能来自随机性而不是参数本身。

%% NSGA-II 主循环 clc; clear; close all; rng(2025); % 固定随机种子,保证结果可复现 % 决策变量边界 [min, max] % x1: 电池容量 (Ah) % x2: 发电机额定功率 (kW) % x3: 电机额定功率 (kW) lb = [20, 2, 5]; ub = [100, 10, 20]; nPop = 100; % 种群规模 nGen = 200; % 迭代代数 nVar = 3; % 决策变量个数 pc = 0.9; % 交叉概率 pm = 0.1; % 变异概率 eta_c = 15; % 交叉分布指数 eta_m = 20; % 变异分布指数 % 初始化种群 population = zeros(nPop, nVar); for i = 1:nPop population(i, :) = lb + rand(1, nVar) .* (ub - lb); end % 计算初始种群目标函数 [obj1, obj2, obj3] = evaluate_HEPS(population); objs = [obj1, obj2, obj3]; for gen = 1:nGen [frontRank, crowding] = non_dominated_sort(objs); selected = tournament_selection(population, frontRank, crowding); offspring = sbx_crossover(selected, lb, ub, pc, eta_c); offspring = polynomial_mutation(offspring, lb, ub, pm, eta_m); [off_obj1, off_obj2, off_obj3] = evaluate_HEPS(offspring); off_objs = [off_obj1, off_obj2, off_obj3]; comb_pop = [population; offspring]; comb_objs = [objs; off_objs]; [population, objs] = environmental_selection(comb_pop, comb_objs, nPop); end

决策变量边界设置也要有依据。电池容量20到100安时对应续航和重量权衡的空间,发电机2到10千瓦覆盖了从纯电池应急到油电并用的不同策略,电机5到20千瓦保证爬升功率需求能被覆盖。边界设太宽会浪费搜索能力,太窄又可能漏掉最优区域,需要根据任务剖面反复试。

3.2 非支配排序与拥挤度距离

非支配排序的作用是把种群分层:第一层是非支配解,第二层是被第一层支配的解,以此类推。这里我使用了标准的分层算法,对每一对个体判断支配关系。拥挤度距离则在同一前沿内计算,目标维度上每个解前后邻居的距离之和,边界解直接设为无穷大保证不会被淘汰。

function [frontRank, crowding] = non_dominated_sort(objs) nPop = size(objs, 1); frontRank = zeros(nPop, 1); crowding = zeros(nPop, 1); dominates = false(nPop, nPop); for i = 1:nPop for j = 1:nPop if i ~= j && all(objs(i, :) <= objs(j, :)) && any(objs(i, :) < objs(j, :)) dominates(i, j) = true; end end end remaining = 1:nPop; currentFront = 1; while ~isempty(remaining) isNonDominated = true(length(remaining), 1); for idx = 1:length(remaining) i = remaining(idx); for jdx = 1:length(remaining) j = remaining(jdx); if i ~= j && dominates(j, i) isNonDominated(idx) = false; break; end end end front_indices = remaining(isNonDominated); frontRank(front_indices) = currentFront; crowd = crowding_distance(objs(front_indices, :)); crowding(front_indices) = crowd; remaining(isNonDominated) = []; currentFront = currentFront + 1; end end

一个常见错误是跨前沿计算拥挤度距离,这样排序意义完全错乱。NSGA-II里的拥挤度比较只发生在同一前沿的个体之间,跨前沿的个体只用层级比较。还有一点,如果某个目标的所有值几乎相同,归一化分母接近零,我在代码里用max - min < 1e-10做了保护,避免除零异常。

3.3 选择交叉变异算子的工程调优

锦标赛选择用的是二元锦标赛:随机抽两个个体,前沿编号小的胜出,前沿编号一样时拥挤度距离大的胜出。这样既保证了收敛压力,又保住了解的分布性。交叉我用的是模拟二进制交叉SBX,变异用多项式变异。这两个都是实数编码的标准算子,代码可以直接复用,但参数需要根据问题调整。

交叉分布指数eta_c直接影响子代与父代的接近程度。eta_c越大,子代越接近父代,种群多样性越差。我一开始设了50,结果帕累托前沿上的点挤成一团,分布非常差,后来调到15才正常。建议eta_c在5到20之间调试,不要上来就取大值。变异概率0.1在3个决策变量时差不多每次迭代有约30%的个体发生变异,多样性基础够用。

环境选择的精英保留策略是NSGA-II的精华:父子种群合并后,先按非支配层级排序,层级小的优先进入下一代,同一层级内拥挤度大的优先。这样保证优秀个体不会在迭代中丢失。

4. 完整目标函数实现与运行结果解读

目标函数是整个仿真框架里最重的部分,它要把飞行任务剖面、部件模型、能量管理策略全部串起来,输出起飞总质量、任务总油耗、电池寿命衰减系数三个目标值。

4.1 目标函数代码逐段解析

function [Mass, FuelCons, BatteryDeg] = evaluate_HEPS(population) n = size(population, 1); Mass = zeros(n, 1); FuelCons = zeros(n, 1); BatteryDeg = zeros(n, 1); g = 9.81; % 重力加速度 rho = 1.225; % 空气密度 S = 1.2; % 机翼面积 (m^2) CD0 = 0.03; % 零升阻力系数 K = 0.05; % 诱导阻力因子 eta_prop = 0.8; % 螺旋桨效率 eta_motor = 0.9; % 电机效率 eta_esc = 0.95; % 电子调速器效率 eta_gen = 0.85; % 发电机效率 LHV = 43.0; % 燃油低热值 (MJ/kg) SFC = 0.35; % 燃油消耗率 (kg/kWh) V_cruise = 30; % 巡航速度 (m/s) T_cruise = 3600; % 巡航时间 (s) T_climb = 300; % 爬升时间 (s) for i = 1:n x = population(i, :); C_bat = x(1); % 电池容量 (Ah) P_gen = x(2) * 1000; % 发电机额定功率 (W) P_motor = x(3) * 1000;% 电机额定功率 (W) V_bat = 50; E_bat = C_bat * V_bat / 1000; % 电池能量 kWh m_bat = E_bat / 0.2; % 电池重量,按200Wh/kg m_gen = P_gen / 1000 / 0.25; % 发电机重量,按0.25kW/kg m_motor = P_motor / 1000 / 1.0;% 电机重量,按1.0kW/kg m_struct = 8.0; % 结构重量 m_payload = 2.0; % 载荷重量 % 能量管理策略 r_climb = 0.8; r_cruise = 0.3; % 爬升功率需求计算 ROC = 3.0; % 爬升率 (m/s) W_guess = (m_bat + m_gen + m_motor + m_struct + m_payload) * g; L_climb = W_guess; D_climb = 0.5 * rho * V_cruise^2 * S * (CD0 + K * (L_climb / (0.5 * rho * V_cruise^2 * S))^2); T_climb_req = D_climb + W_guess * ROC / V_cruise; P_climb = T_climb_req * V_cruise / (eta_prop * eta_motor * eta_esc); % 巡航功率需求 W_cruise = W_guess; D_cruise = 0.5 * rho * V_cruise^2 * S * (CD0 + K * (W_cruise / (0.5 * rho * V_cruise^2 * S))^2); T_cruise_req = D_cruise; P_cruise = T_cruise_req * V_cruise / (eta_prop * eta_motor * eta_esc); % 电池功率与能量消耗 P_bat_climb = r_climb * P_climb; P_bat_cruise = r_cruise * P_cruise; E_bat_total = (P_bat_climb * T_climb + P_bat_cruise * T_cruise) / 1000 / 3600; % 发电机功率与油耗 P_gen_climb = max(0, P_climb - P_bat_climb); P_gen_cruise = max(0, P_cruise - P_bat_cruise); E_gen_total = (P_gen_climb * T_climb + P_gen_cruise * T_cruise) / 1000 / 3600; FuelCons(i) = E_gen_total * SFC; % 电池寿命衰减简化模型 DOD_climb = P_bat_climb * T_climb / (E_bat * 3600); DOD_cruise = P_bat_cruise * T_cruise / (E_bat * 3600); avg_C_rate = max(P_bat_climb, P_bat_cruise) / (E_bat * 1000 / V_bat); BatteryDeg(i) = 0.1 * (DOD_climb + DOD_cruise) * (1 + avg_C_rate / 5); m_fuel = FuelCons(i); Mass(i) = m_bat + m_gen + m_motor + m_struct + m_payload + m_fuel; end end

这里我采用的重量估算公式都是工程经验公式,电池按200Wh/kg能量密度折算,发电机按0.25kW/kg功率密度折算,电机按1kW/kg折算。不同技术水平的部件参数变化很大,实际做方案时需要用真实供应商数据替换。

4.2 从帕累托前沿到工程决策

运行结束后,三维散点图会展示帕累托前沿。通常你会看到一个弯曲的曲面,三个目标之间存在清晰的权衡关系。电池容量大的配置起飞质量必然偏高,但巡航阶段电池承担比例高,油耗和电池寿命衰减相对更低。发电机额定功率大的配置,油耗会上升,但电池衰减会减轻,因为发电机承担了更多基础负荷。

从帕累托前沿上挑点是有方法论的。如果特别看重航程和油耗,应选择目标2最小的点。如果看重寿命和可靠性,选目标3最小的点。如果追求综合均衡,我常用加权距离法选“膝点解”,这个点一般是曲率最大的位置,代表整体性能最优的折中方案。

weights = [0.4, 0.3, 0.3]; normalized = (pareto_objs - min(pareto_objs)) ./ (max(pareto_objs) - min(pareto_objs)); dist = normalized * weights'; [~, best_idx] = min(dist); best_solution = pareto_solutions(best_idx, :); best_obj = pareto_objs(best_idx, :);

这里的权重是决策者偏好,是在帕累托前沿解集里选点,和直接把三个目标加权合并成单目标优化有本质区别。后者从搜索一开始就限定了可能性空间,前者则保留了完整的设计权衡信息。

5. 常见问题排查与工程避坑实录

我把这个项目从零跑到最终结果的过程中踩了不少坑,挑几个典型的列出来,省得后来者重复走弯路。

5.1 结果全部聚集在边界上

如果优化完发现所有解都挤在决策变量边界附近,先别急着调算法。按这个顺序排查:检查目标函数在边界附近是否有突变或数值溢出,交叉验证某个目标是否完全不受某些决策变量影响,最后检查约束处理。NSGA-II本身不强制约束,如果目标函数里没有做任何约束惩罚,算法会把明显不可行的配置也当成好解。工程上的约束比如功率需小于发电机最大功率、电池能量必须覆盖任务剖面等,都应该在目标函数里加惩罚项,或者直接让超出约束的解产生极大的目标值,把搜索引回可行域。

5.2 帕累托前沿点太少且分布不均

多半是拥挤度距离没算对,或者选择压力过大导致种群多样性丧失。拥挤度距离必须同一前沿内计算,边界个体设为大值。选择压力方面,增大变异概率和降低交叉分布指数可以提升多样性,但要小心破坏收敛性。我从eta_c=15调eta_c=5之后,前沿分布明显更均匀,但收敛代数也变长了,需要根据可接受的计算时间做平衡。

5.3 每次运行结果差异巨大

多目标进化算法本质有随机性,这不一定代表代码有Bug。但如果你需要可复现的实验结果,固定随机种子是硬性要求。rng(2025)这行代码让每次运行结果一致,调参时可以精确对比。不同随机种子之间的结果差异本身也能反映算法的稳定性,如果换了种子前沿形状变化很大,说明种群规模或者迭代代数不够,需要调大。

5.4 巡航功率算出来明显不合理

我在调试时单独写了一个测试脚本,随机生成一组决策变量,调用目标函数后打印所有中间变量:电池重量、发电机重量、各阶段功率、燃油重量、总重量。如果巡航功率只有几十瓦或者燃油重量出现负数,那一定是模型有Bug,而不是算法的问题。有个非常隐蔽的错误:螺旋桨效率和电机效率连乘时忘记除以,我一度把功率需求算得偏小,结果整个质量估算偏低,优化出的配置全都不满足真实任务需求。打印中间变量逐个对照物理常识是排查这类问题最高效的方式。

5.5 目标计算性能瓶颈

如果以后把种群规模扩展到500甚至1000,或者任务剖面细化成十几个阶段,for循环逐个个体算目标函数会明显变慢。建议把目标函数向量化,一次处理整个种群,虽然可读性稍微差一点,但性能提升显著。以我这个项目的规模看,100个个体跑200代,几分钟就能完成,不需要过早优化性能。

6. 向其他混动系统迁移这套框架

这套NSGA-II加Matlab的框架可以很自然地迁移到其他混合动力场景,比如电动垂直起降飞行器、混合动力汽车、甚至船舶混合动力系统。迁移的核心不是算法代码,而是系统模型和目标函数。

步骤很清晰:首先替换推进架构模型,把电池、发电机、电机之间的能量流关系按新系统的拓扑结构改掉。然后调整目标函数,匹配新系统的关键性能指标,比如汽车可以加百公里加速时间目标,船舶可以加排放目标。接着重新设置决策变量边界,根据实际可用部件供应商的规格范围来确定搜索空间。最后重新标定重量和效率参数,这是最花时间的环节,每个部件的技术参数都必须基于可靠数据。

算法层面的主循环、非支配排序、拥挤度计算、锦标赛选择、SBX交叉、多项式变异这些模块完全不用改。这也从侧面反映了NSGA-II这类经典多目标算法的通用性。

我个人更建议不要一上来就追求高保真模型。第一次跑通整个框架时,用简化模型完全足够,先把算法链路走通、确认帕累托前沿形状合理,再去逐步增加复杂度,比如加入电池热模型、变高度变速度飞行剖面、螺旋桨变距控制等。渐进式搭建比一次堆出巨型模型然后陷入调试泥潭高效得多。我在这套模拟项目X上就是先跑通简化模型,再逐步完善,整个开发过程节奏比较顺畅。

最后再强调一句我在实际调试中不断验证的体会:多目标优化的最大价值不在于算法本身跑得多漂亮,而在于目标函数里是否真实还原了物理系统的耦合关系。当每个帕累托前沿上的点都对应一组有明确工程含义的设计参数时,这份优化结果才能真正支撑设计决策。如果你正准备把NSGA-II用到类似的飞行器能源系统优化中,这套代码可以作为一个可靠起点,在此基础上针对你的机型改参数、改剖面、改约束,就能很快得到有意义的结果。

RELATED READING

延伸阅读

更多一线实战笔记与深度复盘,助您持续精进