前阵子帮朋友调试一个调度优化项目一座中型水电站和一片光伏电站打捆外送需要排日前24小时的发电计划。他原来的方案是光伏满发、水电兜底看起来很简单但实际曲线很难看——白天光伏出力冲高水电只能压着傍晚光伏骤降水电要在一两个小时内补上出力缺口送出功率的爬坡率远超电网允许范围。更麻烦的是这个调度要同时照顾发电收益、弃光电量、出力平稳性三个指标改任何一个时段都可能牵动全局。显然这不是一个能靠经验手调的优化问题而是典型的多目标优化问题。我最后选定的是非支配排序遗传算法NSGA-II求解多目标水光互补优化调度模型全部用Python实现。这份记录会把从问题建模、算法原理到代码落地、调参踩坑的完整过程写清楚适合正在做电力优化调度、水库运行规划或者刚接触多目标优化想找个实际工程案例练手的朋友参考。1. 问题拆解水光互补调度到底在优化什么1.1 为什么要把水电和光伏捆在一起调度水电的最大优势是调节快水轮机组从开机到满发通常几分钟内完成库容可以起到储能作用白天少发存水晚上多发补位。光伏则正好反过来白天出力大、夜间为零而且受云层影响时刻波动单独并网时电网要时刻准备旋转备用。两者打捆后可以让水电去追踪光伏的波动用水电的确定性平抑光伏的不确定性同时共享送出通道提高线路利用率。这个逻辑说起来简单落到数学上就很复杂水电到底在哪几个小时发、发多少既要考虑来水能不能支撑又要考虑库容会不会越界还要考虑水头变化对出力的影响。一个时段的决策会通过水量平衡传导到后面所有时段所以它不是一组独立的24个变量而是带时间耦合的约束优化问题。1.2 三个目标函数收益、弃光、平稳性我在这个项目里选了三个典型目标全部按最小化处理发电收益负值 f1 -∑ price_t × (Ph_t Ppv_t)价格取分时电价。加负号是为了统一成最小化算法内部不需要区分方向。弃光电量 f2 ∑ (Ppv_ava_t - Ppv_t)其中 Ppv_ava_t 是t时段光伏可用出力Ppv_t 是实际消纳出力只能小于等于可用值。送出功率波动 f3 ∑ |P_out_t - P_out_t-1|其中 P_out_t Ph_t Ppv_t。这个目标直接对应电网对爬坡率的限制。这三个目标天然冲突想让弃光少就要多消纳光伏但光伏集中在白天这意味着水电最好在白天少发、把库容留给光伏可这样傍晚补位时段的出力爬坡就会很猛想让曲线平稳则要牺牲一部分高电价时段的出力。多目标优化的价值就在这儿不是硬拍一个权重而是把这一组取舍关系完整算出来让决策的人事后按偏好挑方案。工程上还会加第四个目标比如弃水量最小、调度期发电量最大、水位偏差最小。目标不是越多越好3到4个比较合适目标太多NSGA-II会因为帕累托前沿维度过高而难以均匀采样这时候得考虑NSGA-III那类处理高维目标的算法。1.3 约束条件要列全否则结果没法用优化模型的约束我按硬约束优先、软约束惩罚的思路处理核心约束有四类。第一是水量平衡。水库蓄水量 V_t 的递推关系是 V_t1 V_t I_t - Qh_t - Qspill_tI_t 是入库流量Qspill_t 是弃水。这里注意弃水不是决策变量而是计算结果——如果某时段即使发电流量拉到上限库容还是会超过最大库容 V_max多出来的水只能弃掉。第二是边界约束。库容不能低于死库容 V_min也不能超过正常蓄水位对应的 V_max发电流量 Qh_t 介于 Qh_min 和 Qh_max 之间水电出力 Ph_t 介于技术出力下限和装机上限之间。水电出力我用的简化公式 Ph_t K × Qh_t × H_tK是综合效率系数H_t是净水头。净水头不能当成常数这个坑后面我会专门讲。第三是末端水位约束。调度期末库容 V_T 要大于等于某个目标值比如维持正常蓄水位的80%。没有这个约束算法会把水在最后几天全发完第二天就没法接续了这在滚动调度里是必须考虑的。第四是光伏出力约束。实际调度中如果要弃光Ppv_t 可以取小于可用值的数如果不允许弃光那就直接令 Ppv_t Ppv_ava_t把弃光目标从模型里拿掉。我这个算例里把弃光当成目标本身所以人为允许调度去选择弃掉一部分光伏来换取更平稳的外送曲线。还有一个隐含约束是 P_out_t 要满足电网允许的最大外送功率和爬坡率限值。爬坡率限值如果写成 |P_out_t - P_out_t-1| ≤ 斜坡上限它其实和波动目标 f3 是同一件事区别只是一个是硬约束一个是目标。我这里选择把它放在目标函数里让算法自己权衡爬坡和经济性而不是一刀切限制死。2. 为什么选NSGA-II多目标优化算法选型逻辑2.1 先搞清支配和帕累托前沿多目标优化里没有全局最优解这个概念只有非支配解。还是以最小化为例解A支配解B的条件是A在所有目标上都不比B差且至少有一个目标严格优于B。一个种群中所有不被任何解支配的个体构成第一前沿也就是帕累托前沿被第一前沿支配、但不被其它解支配的个体构成第二前沿依此类推。帕累托前沿上的解没有谁绝对优于谁它们只是在目标空间里占据了不同的取舍位置。放在我的调度问题里前沿上的一个解可能是收益最高但弃光偏多、波动较大另一个可能是曲线最平稳但总收益少了一截。调度的任务就是先把这一整条前沿找出来再按实际偏好挑一个点。2.2 NSGA-II的三个关键机制拆解Deb等人在2002年提出的NSGA-II核心机制就是三件事。第一是快速非支配排序。一个种群有N个个体两两比较目标向量得到每个个体被谁支配、支配谁的信息然后按支配关系把种群分层。这里有个重要细节排序只看目标值不看约束是否满足。约束处理要靠目标函数里的惩罚项或单独修复排序环节不需要考虑可行性否则很容易丢掉有潜力的边界个体。第二是拥挤度距离。同一前沿层内的个体按每个目标维度分别排序相邻个体的目标差除以该目标取值范围归一化后累加得到一个周围有多空的指标。端点个体距离直接设为无穷大保证前沿的两个端点不被淘汰这样前沿在目标空间里能展开到足够宽。第三是精英保留。每一代把父代和子代合并成2N规模先按前沿层数从小到大排层数小的优先入选如果某一层放不下全部个体就用拥挤度排序保留距离大的个体。这样上一代表现最好的解不会因为随机交叉变异被冲掉收敛性比老版NSGA强很多。2.3 为什么不选加权法、MOPSO或NSGA-III做这个项目之前我把主流的多目标方法都过了一遍。加权法最省事但有两个硬伤权重本质上是在拍脑袋决定目标之间的相对重要性而且如果帕累托前沿非凸加权法无论怎么扫权重都找不到前沿凹进去的部分。这对非线性的水光调度模型来说是个致命伤。MOPSO多目标粒子群在连续变量问题上通常收敛很快但粒子群多了惯性权重、个体学习因子、群体学习因子、速度上限四个额外参数每个参数对结果都敏感。我的经验是工程问题上粒子群调参的时间成本比遗传算法高不少而且对约束问题的处理没有遗传算法顺手。NSGA-III是处理4个以上目标的高维版本水光调度一般3个目标就够用不上。综合下来选了NSGA-II实现难度中等、鲁棒性好、社区资料多出了问题也容易排查。3. Python实现从目标函数到算法主循环3.1 编码设计与种群初始化决策变量我定义为两组合并成一个长度为48的实数向量前24位是归一化的发电流量比例范围[0,1]实际流量再乘一个最大流量系数后24位是各时段的光伏弃光比例 alpha_t范围也是[0,1]实际消纳光伏等于 Ppv_ava_t × (1 - alpha_t)。为什么用比例而不是直接用实际值因为交叉和变异算子在[0,1]区间里操作最稳不会产生非法流量值和负光伏出力也方便后续的修复操作。种群用NumPy存成 (N, 48) 的二维数组一行就是一个个体。初始化不能全用均匀随机数。我试过均匀随机生成的个体大部分末端库容严重越界种群起点就很差。改进方法很朴素生成随机数后做一次流量平移修复把24时段的平均流量往满足末端库容的方向调整让初始个体至少是边界可行的。这一小步能让后续进化快很多推荐直接用。3.2 目标函数和约束评估怎么写评估函数是每次迭代调用最频繁的部分要把水量平衡、水头、出力的计算写成一个清晰的函数。下面是我用的评估框架def evaluate(pop, inflow, pv_ava, price, params): N, D pop.shape T len(pv_ava) Qh pop[:, :T] * params[Qh_max] # 实际发电流量 alpha pop[:, T:] # 弃光比例 0~1 Ppv (1 - alpha) * pv_ava # 实际消纳的光伏出力 V np.zeros((N, T 1)) V[:, 0] params[V0] Ph np.zeros_like(Qh) for t in range(T): Z params[a] * V[:, t] params[b] # 水位 H Z - params[tail_level] # 净水头 Ph[:, t] params[K] * Qh[:, t] * H # 水电出力 V[:, t1] V[:, t] inflow[t] - Qh[:, t] spill np.maximum(0, V[:, t1] - params[V_max]) V[:, t1] - spill # 超出库容部分弃水 P_out Ph Ppv f1 -np.sum(price * P_out, axis1) f2 np.sum(alpha * pv_ava, axis1) # 弃光电量 diff P_out[:, 1:] - P_out[:, :-1] f3 np.sum(np.abs(diff), axis1) # 惩罚项末端库容不足、运行期库容越界 penalty 1e6 * (np.maximum(0, params[V_end_req] - V[:, T]) np.maximum(0, params[V_min] - V.min(axis1)) np.maximum(0, V.max(axis1) - params[V_max])) return f1 penalty, f2 penalty, f3 penalty这段代码有几个可以注意的点。一是惩罚项加在三个目标上而不是单独记一个约束值这样非支配排序时可行解和不可行解之间的比较会自然把不可行解往后压但又不会完全丢掉那些目标非常优秀但轻微越界的个体给修复留了余地。二是弃光通过 alpha 变量进入决策模型才能真正在多喘光和少弃光之间做取舍如果直接把Ppv固定为可用值f2恒等于0前沿会退化整个多目标就名存实亡了。实际项目中光伏要不要弃通常由调度指令决定。如果上游已经下达了光伏优先消纳的指令那模型里就不该给弃光自由直接固定消纳率如果允许弃光调峰才用上面的 alpha 写法。3.3 快速非支配排序核心代码排序直接照搬Deb论文思路用两层循环比较支配关系。虽然复杂度是O(MN²)M是目标数N是种群数但N一般只有100到200完全够用。import numpy as np def dominates(a, b): # 所有目标均为最小化 return all(x y for x, y in zip(a, b)) and any(x y for x, y in zip(a, b)) def fast_non_dominated_sort(obj_matrix): N obj_matrix.shape[0] S [[] for _ in range(N)] # 被个体i支配的个体集合 n np.zeros(N, dtypeint) # 支配个体i的个体数量 fronts [[]] for i in range(N): for j in range(N): if i j: continue if dominates(obj_matrix[i], obj_matrix[j]): S[i].append(j) elif dominates(obj_matrix[j], obj_matrix[i]): n[i] 1 if n[i] 0: fronts[0].append(i) k 0 while fronts[k]: q [] for i in fronts[k]: for j in S[i]: n[j] - 1 if n[j] 0: q.append(j) k 1 fronts.append(q) return fronts[:-1]一个容易写错的地方是 S 和 n 的更新顺序必须先收集完第一前沿再逐层处理不能在外层循环里边算边更新否则层数会乱。3.4 拥挤度距离维持前沿多样性的关键同一层内的个体按每个目标维度分别排序相邻个体的距离差除以该目标取值范围累加得到拥挤度def crowding_distance(front_idx, obj_matrix): dist {idx: 0.0 for idx in front_idx} m obj_matrix.shape[1] for obj in range(m): order sorted(front_idx, keylambda idx: obj_matrix[idx, obj]) dist[order[0]] np.inf dist[order[-1]] np.inf fmin obj_matrix[order[0], obj] fmax obj_matrix[order[-1], obj] if fmax - fmin 1e-12: continue for j in range(1, len(order) - 1): dist[order[j]] (obj_matrix[order[j1], obj] - obj_matrix[order[j-1], obj]) / (fmax - fmin) return dist注意阈值 1e-12如果某个目标在所有个体上取值都一样分母为0必须先跳过。我在初版没加这个判断结果算出一大片无穷大选择操作直接失效。3.5 SBX交叉与多项式变异实数编码对应的经典算子是模拟二进制交叉SBX和多项式变异。SBX的思想是让子代围绕父代产生同时通过分布指数 eta_c 控制子代贴近父代的程度eta_c 越大子代离父代越近。def sbx_crossover(p1, p2, eta_c15, prob0.9): c1, c2 p1.copy(), p2.copy() mask np.random.random(p1.shape) prob u np.random.random(p1.shape) beta np.where(u 0.5, np.power(2*u, 1/(eta_c1)), np.power(2*(1-u), -1/(eta_c1))) c1[mask] 0.5 * ((1beta)*p1 (1-beta)*p2)[mask] c2[mask] 0.5 * ((1-beta)*p1 (1beta)*p2)[mask] return np.clip(c1, 0, 1), np.clip(c2, 0, 1) def polynomial_mutation(child, eta_m20, prob0.1): r np.random.random(child.shape) delta np.where(r 0.5, np.power(2*r, 1/(eta_m1)) - 1, 1 - np.power(2*(1-r), 1/(eta_m1))) child child (r prob) * delta return np.clip(child, 0, 1)变异概率指每个维度的变异概率48个维度每个维度独立判断。不要把整条染色体的变异概率设成0.1那等于整条曲线只有四五个点会动前期种群多样性会非常差。3.6 精英选择与主循环主循环按NSGA-II标准流程走随机初始化并修复种群 P评估目标值二元锦标赛选择父代比较规则是先比前沿层层数小的赢层数相同比拥挤度距离大的赢对选出的父代做SBX交叉和多项式变异生成子代 Q合并 R P ∪ Q规模2N评估目标值对 R 做快速非支配排序从第一前沿开始依次填充下一代填到某一层时如果剩余名额不足用拥挤度排序保留距离大的个体重复2到6直到达到最大代数。刚接触NSGA-II的朋友问我怎么检查实现对不对我的建议是先拿ZDT1、ZDT2这种标准测试函数跑一遍比对已知前沿形态。能通过标准测试再上水光调度模型能省一半的调试时间。标准测试就是验算台代码写出来是一回事代码写对是另一回事。4. 仿真结果分析帕累托前沿长什么样4.1 算例参数与数据准备我做的一个典型算例参数如下参数取值调度周期24小时步长1小时水电装机30 MW光伏装机50 MW最大发电流量45 m³/s死库容 / 正常库容80 / 200 万m³初始库容150 万m³末端库容要求≥ 120 万m³入库流量枯水期典型日径流取10~20 m³/s光伏可用出力夏季晴天带云曲线中午接近45 MW光伏可用出力的数据如果是实际项目最好用当地气象源的历史数据或者调度部门给的预测值没有数据的话自己生成一条平滑偏态曲线也够验证算法逻辑。入库流量同理短期预报数据更可靠。价格取分时电价峰段(17:00-21:00) 1.0元/kWh平段(08:00-17:00) 0.6元/kWh谷段(21:00-08:00) 0.3元/kWh。水电和光伏的度电成本有差异但算例里我统一用上网口径比较发电端收益结构不细分电源成本。4.2 帕累托前沿的典型形态种群规模100、进化300代跑完把三个目标两两画散点图能看到一个非常典型的弯曲前沿面。收益和弃光之间的趋势是收益越高的解弃光越少因为多发的电都是有价值的但同时收益越高的解出力波动往往越大因为高收益要求晚高峰满出力、白天又要给光伏让路两部分衔接处的爬坡自然就猛。如果三目标前沿在3D图里看应该是一张略有弧度的弯面而不是一条直线。如果退化成一条直线说明某个目标基本是其它两个目标的线性组合如果缩成一个点说明种群多样性坏了。我一般还会算一个超体积指标HVhypervolume定量评价HV随时间上升的曲线是判断收敛最直观的工具。4.3 结果的合理性验证算出来的方案合不合理不能只看目标值。我习惯再做一个人工规则方案做对比规则方案就是光伏满发水电在晚高峰前蓄水、在晚高峰满发补位同样计算三个目标值然后把这个方案的目标点画到帕累托前沿图里。正常情况下规则方案的点应该落在前沿下方或里面也就是说NSGA-II找出的非支配解至少在三个目标上都优于或持平于人工规则。我这次跑出来的结果里规则方案明显被前沿支配同样波动水平下NSGA-II解的收益高出约6%同样收益水平下波动小了一半。这说明模型和算法没有空转是真的在目标空间里掘出了比经验更优的调度方式。4.4 从解集里挑最终调度方案帕累托前沿只是一堆候选解最终要用一个。我推荐先用熵权TOPSIS从目标矩阵里算综合排序选出折中点如果现场有明确偏好也可以在前沿上直接加筛选条件。这个事后加偏好的能力正是用帕累托方法比加权法舒服的地方上午拿到的前沿下午想改成弃光率小于5%且波动尽量小直接筛不用重新跑算法。筛选出方案后我会把解对应的24小时库容曲线、水电出力曲线、光伏实发曲线和送出功率曲线全部画出来确认没有明显的时序异常。特别是库容曲线它必须是平滑的、不越界的如果出现锯齿多半是惩罚项还是偏软需要调大惩罚系数或者改用更严格的修复策略。5. 踩坑记录与调参经验5.1 惩罚函数系数是最大的坑初版实现里我把库容越界惩罚设成和发电收益同一个量级结果种群几乎全是不可行解前沿根本展不开。后来改成了动态惩罚记录当前代的平均约束违反量惩罚系数随代数缓慢放大前期允许个体试探性地越界后期严格压制。效果立竿见影。更推荐的做法是修复式初始化在种群生成后把不满足末端库容的个体按比例整体调整发电流量让它刚好落在要求附近。这样大部分个体先天就是可行的惩罚项只处理运行时段的短期越界种群可行性会高很多。5.2 水头不能当常数处理这个坑我印象很深。初版把净水头固定成常数算法给出的最优调度会把水库在中后期放得很低因为算法以为水头始终不变想尽量在高水位时段多拉流量。但一旦考虑水头随库容下降中后期水头低、同样的流量发电量会缩水原来那个方案就不再最优。加上水位-库容线性近似之后调度曲线明显更均衡结果也更符合常识。工程上更准确的做法是直接查水位-库容曲线和水头-尾水位关系表但这需要厂家数据或调度台账。自己搭模型阶段用 a×Vb 的线性近似完全够用关键是水头会随库容变化这个机制必须在模型里体现否则解的物理意义会失真。5.3 收敛判断别只看代数NSGA-II没有单一收敛指标。我自己看两个东西一是HV指标的增量曲线连续50代变化小于阈值就算收敛二是看目标空间里前沿的展开面积如果前沿还在明显拉长说明还有潜力继续跑。200到500代是水光调度模型的常见区间但具体代数取决于问题规模和约束松紧不能照搬别人的数值。种群规模方面我建议先试100如果前沿空档明显两个解之间距离很大就加到150或200。小种群跑得快但前沿稀疏后处理选方案时会很尴尬大种群跑得慢但前沿丰满选方案方便。两三个目标、48个决策变量的规模200×300代在普通笔记本上也就几分钟完全跑得起。5.4 随机种子和重复实验多目标遗传算法本质是随机算法单次运行的结果不能作为结论。我在项目里每个参数组合跑5次把5次的帕累托前沿合并后再做一次非支配排序取合并后的第一前沿作为最终前沿。这样可以抹平单次随机波动也让最终录取的方案更有说服力。写报告或论文时把5次运行的目标值区间列出来比只画一个最佳前沿可信得多。另一个我自己常用的检查是把最优个体的决策变量还原成调度曲线请熟悉现场运行的同事看一眼看水位过程线是否符合实际运行规则。算法结果合理但不符合老师傅直觉时通常不是算法错了而是模型漏了某条约束——比如下游生态流量下限、机组最小开停机时间。把这条约束补进去再跑结果很快会变得像样。跑这个项目最大的体会是多目标调度模型的难点一直不在算法本身而在目标函数和约束的物理还原度。NSGA-II的代码结构就那些但水量平衡、水头变化、末端水位这些机制每补一个解的质量就会上一个台阶。如果你也要做类似的水光互补或水库优化调度优先把时间花在模型细节上算法版本反而是最不用担心的一环。