风光柴储微网多目标优化调度：NSGA-II与MATLAB实现详解

做电力系统优化调度的人应该都有共识：风光柴储微网的多目标优化调度，最难的不是模型本身，而是“怎么把模型写进代码、跑出结果、再让结果经得起推敲”。我之前被这个题目折腾了挺久，一开始拿网上的通用框架跑，总是要么成本降不下来，要么储能SOC直接越界，后来干脆自己照着原理从零写了一套MATLAB代码，用多目标优化算法（NSGA-II）去求解24小时机组出力计划，才算把经济性和环保性这两个矛盾目标理顺。这篇文章就是我当时的那套思路和代码结构，从数学模型搭建、算法选型，到MATLAB实现细节、仿真结果解读，再到踩过的坑，一次性讲完。

所涉及的问题适合正在做微网调度、综合能源系统优化或者正在写多目标优化方向论文的同学参考。我会尽量把每个“为什么这么写”都说清楚，因为光贴代码没意义，能照着改、照着调，才是这套东西的价值。

1. 问题建模：风光柴储微网调度的数学模型怎么搭

1.1 系统结构与目标函数拆解

先明确研究对象。所谓“风光柴储微网”，系统内部最核心的单元就是四类：风力发电、光伏发电、柴油发电机、储能电池。微网与大电网之间还留了一个联络线，可以买电也可以卖电。调度周期取典型日24小时，时间步长1小时，决策结果就是未来24小时内各个可控单元的出力计划。

在这个系统里，风电和光伏属于不可控电源，它们的出力上限由预测曲线决定，调度模型不能要求它们超出预测值出力，但可以选择弃风、弃光，也就是实际出力可以小于预测出力。柴油机和储能是可控单元，柴油机可以调节出力，储能可以充放电。上级电网交换功率也可以作为调节手段，但通常受到联络线容量约束。

目标函数一开始我列了三个：运行经济性、污染物排放量、以及电能质量相关的指标。但为了不让问题过于复杂，也为了能画出一张清晰的帕累托前沿图，我最终保留了前两个目标，第三个目标放在约束条件里间接考虑。

第一个目标是综合运行成本，公式上一般拆成四块：

• 柴油机燃料成本：采用二次耗量特性曲线，即 a_i + b_i * P_i(t) + c_i * P_i(t)^2，这是经典做法；
• 设备运行维护成本：风电、光伏、柴油机、储能各自乘以单位运维系数；
• 柴油机启停成本：因为模型里用连续出力变量，启停成本可以用一个简化表达式近似处理，比如每次状态从0到非0时增加一个固定费用；
• 与电网的购售电费用：购电价高于售电价，所以目标函数鼓励自发自用。

第二个目标是污染物排放量，主要是柴油机燃烧产生的CO₂、SO₂、NOₓ，统一折算成等效排放量，然后乘以排放系数。风电、光伏、储能在运行阶段基本零排放，所以排放目标实际上只跟柴油机和电网购电的间接排放有关。这样做虽然简化了，但已经能反映出“多目标”的本质：成本最小化通常会让柴油机多出力，而排放最小化则希望压低柴油机出力、多用新能源和储能，两个目标方向相反，这就构成了非劣解集。

1.2 约束条件：功率平衡、机组出力、储能SOC都不能踩线

约束条件才是调度模型真正花时间的地方。我归纳了一下，至少包含以下五类：

1. 功率平衡约束。任意时刻 t，光伏出力加风电出力加柴油机出力加储能放电功率加购电功率，必须等于负荷加储能充电功率加售电功率。这是等式约束，也是整个模型最核心的一条线。实际编程时，我基本不会把这条等式的所有项都作为独立决策变量，那样等式约束很难严格满足，后面会说残差变量的处理技巧。

2. 柴油机出力上下限约束。单台柴油机有最小稳定技术出力，通常取额定功率的30%，也有取20%的，具体看机型。如果调度周期内包含停机状态，还要考虑连续启停时间约束。简化处理时，我会直接要求各时段出力要么为0，要么在上下限区间内。

3. 爬坡约束。柴油机从一个时段到下一个时段的出力变化量不能超过爬坡速率。这个约束容易在代码里漏掉，跑出来的调度曲线经常出现相邻时段从10kW直接跳到50kW的极端情况，现场肯定不允许。

4. 储能SOC约束与充放电功率约束。储能电池的SOC按时间递推：SOC(t) = SOC(t-1) + (η_ch * P_ch(t) - P_dis(t)/η_dis) * Δt / E_cap。SOC有上下限，通常取0.1到0.9，避免过充过放。充放电功率也有最大值限制，并且同一时段不能同时充电和放电。这个约束如果不用逻辑变量，处理起来有点麻烦，我倾向于在个体修复环节强制写入：当天充放电状态由符号判定，只允许其中一个方向。

5. 联络线功率约束。与电网的交换功率不能超过联络线允许容量，一般取微网总负荷峰值的30%到50%，具体看上级配电网给的接入方案。

这些约束里，功率平衡是等式约束，其他都是不等式约束。处理方式有罚函数法、修复法、以及将部分变量作为残差量的方法，我后面会重点讲残差变量的思路，因为这是让NSGA-II快速收敛的关键。

2. 求解思路与算法选型：多目标优化算法怎么选

2.1 为什么选NSGA-II而不是加权求和

一开始很多人会想：把成本和排放各自乘个权重，加起来变成一个单目标，然后就能用现成的线性规划或遗传算法去求解，不是更简单吗？确实简单，但有个本质问题：权重系数怎么定？成本和排放量纲不同，权重本质上代表决策者对两个目标的偏好，这个偏好换一个人可能就完全不一样。更麻烦的是，权重不均匀扫描时很容易漏掉非凸帕累托前沿上的部分解。

所以我直接上多目标进化算法，用的是NSGA-II。选择它的理由很朴素：第一，MATLAB生态里NSGA-II的资料最多，出图方便，改造成本低；第二，它基于非支配排序和拥挤度距离，既保证收敛性又保证多样性，一次运行就能得到一组帕累托前沿解；第三，不需要提前处理量纲差异，哪怕成本是几万元、排放是几百公斤，直接扔给算法就行。

NSGA-II的核心流程是：初始化种群，然后进入快速非支配排序，把种群分成多个前沿等级；同一前沿内部再用拥挤度距离排序，保证解的分散程度；然后通过锦标赛选择、模拟二进制交叉、多项式变异生成子代，父子合并后再排序截断。循环迭代，直到满足终止条件。

我也可以提一下其他算法作为备选：MOPSO（多目标粒子群）收敛速度快，但解的分布有时候不如NSGA-II均匀；SPEA2在多样性维护上很有特色，但实现代码稍复杂。如果你用的是新版MATLAB，Global Optimization Toolbox里其实自带多目标遗传算法solve函数与gamultiobj，不过它封裝程度高，对想深入理解机理的人来说，还是自己写一遍收获更大。

2.2 决策变量编码与残差变量思路

决策变量怎么选，直接决定优化的搜索难度。最容易想到的做法是：把每个时段的柴油机出力、储能充放电功率、电网交换功率都作为变量。比如24小时、4类决策变量，那就是96维搜索空间。这么做也不是不行，但约束处理会很痛苦，功率平衡等式几乎不可能精确满足，罚函数惩罚系数还特别难调。

我采用的思路是“残差变量法”。核心思想是：功率平衡等式里的某一项不参与随机编码，而是在目标函数计算时用其他变量反推出来。具体来说，我把决策变量设计成两组：

• 24小时的柴油机出力序列；
• 24小时的储能充放电功率序列。

然后电网交换功率 P_grid(t) 由功率平衡公式直接计算得出，不再作为独立变量：

P_grid(t) = P_load(t) - (P_pv(t) + P_wind(t) + P_discharge(t) - P_charge(t) + P_dg(t))

这样一来，只要P_grid(t)落在联络线约束范围内，功率平衡等式就自动满足，不需要额外的罚函数项。如果P_grid(t)越界，那么再对该个体的出力做小幅修正，或者直接把这个越界量作为惩罚加进去。这个“残差变量”的做法在微网调度代码里非常实用，强烈推荐。

储能同时充放电的问题，我也在修复阶段处理：如果某个时段P_ch和P_dis都大于0，就按净功率折算成单一方向，并对SOC递推公式重新计算。这样做虽然牺牲了一部分理论精确性，但换来的是约束的严格满足，对进化算法来说更友好。

3. MATLAB代码实现：从零搭一套可跑的调度程序

3.1 数据初始化与主流程框架

MATLAB代码我拆成了几个文件，职责分开，这样调试时不用在一个几千行的脚本里上下翻：

• main_optimize.m：主程序，负责参数初始化、调用优化、输出结果；
• obj_fun.m：目标函数，输入一个个体，输出成本和排放两个目标值；
• constraint_check.m：约束校验，返回越限量；
• nsga2_operator.m：NSGA-II核心算子（非支配排序、拥挤度计算、选择、交叉、变异）。

先看数据初始化部分。所有基础参数我都用结构体para统一管理，避免到处散落全局变量：

% 基础参数 para.dt = 1; % 时间步长 1小时 para.T = 24; % 调度周期 para.Ndg = 2; % 柴油机台数 para.pv = [%; ...]; % 24小时光伏预测出力 (kW) para.wind = [%; ...]; % 24小时风电预测出力 (kW) para.load = [%; ...]; % 24小时负荷预测 (kW) para.Ecap = 100; % 储能容量 kWh para.SOC0 = 0.5; % 初始SOC para.SOC_min = 0.1; para.SOC_max = 0.9; para.Pch_max = 25; % 最大充电功率 kW para.Pdis_max = 25; % 最大放电功率 kW para.Pgrid_max = 80; % 联络线最大交换功率 kW para.DG_min = [15; 15]; % 各柴油机最小出力 para.DG_max = [50; 50]; % 各柴油机最大出力 para.ramp = 20; % 爬坡约束 kW/小时 para.eta_ch = 0.9; % 充电效率 para.eta_dis = 0.9; % 放电效率

光伏、风电、负荷这组数据，如果是做真实项目，应该来自SCADA系统或气象预测服务；如果是复现论文，可以直接取典型日的公开数据，或者自己按日曲线形状造一组。我在案例里用的是某微网示范工程的夏季典型日数据，光伏白天明显有驼峰状，风电夜间出力大，负荷早晚各有一个峰值。

主程序框架相对固定：

numIndividual = 100; % 种群规模 numGen = 200; % 迭代次数 numVar = 48; % 决策变量维度：2台柴机24小时 + 储能24小时 = 72，示例 % 这里实际按动态变量数目设计，见后文说明

我上面只写了48维，但实际代码里决策变量数量需要跟柴油机台数和调度周期严格对应。比如2台柴机24小时出力序列占48维，储能充放电又占48维，总维数是96维。初代种群用随机数在上下限之间均匀生成，然后用“残差变量法”和约束修复逻辑得到合法个体，再进行进化迭代。

3.2 目标函数与约束处理的关键代码

目标函数代码是整套程序的心脏。我简化后的版本长这样：

function [cost, emission] = obj_fun(x, para) % 解码决策变量 numDg = para.Ndg; P_dg = x(1 : numDg*para.T); P_dg = reshape(P_dg, numDg, para.T); P_bat = x(numDg*para.T+1 : end); P_ch = max(0, P_bat); % 充电功率 P_dis = max(0, -P_bat); % 放电功率，负数表示充能 % 功率平衡计算电网交换功率 P_grid = para.load - (sum(P_dg,1) + para.pv + para.wind + P_dis - P_ch); % 成本目标：燃料 + 运维 + 购电 + 启停 % 燃料成本按二次耗量曲线 cost_fuel = 0; for i = 1 : numDg cost_fuel = cost_fuel + sum(para.a(i) + para.b(i).*P_dg(i,:) + para.c(i).*P_dg(i,:).^2); end % 运维成本 cost_om = sum(para.k_pv.*para.pv + para.k_w.*para.wind + ...); % 购电/售电成本 cost_grid = sum(max(P_grid,0).*para.price_g + min(P_grid,0).*para.price_s); % 启停成本 cost_start = para.cost_start .* sum(diff([zeros(1,numDg); P_dg']~=0, 1) == 1); cost = cost_fuel + cost_om + cost_grid + cost_start; % 排放目标：简化等效排放系数 emission = sum(sum(para.emis_para * P_dg)) + sum(max(P_grid,0).*para.emis_grid); end

注意代码里我用了reshape把个体向量还原成矩阵，这是MATLAB里处理多维决策变量最常见的方式。奥妙在于把“柴油机出力”和“充放电功率”分开编码，但充放电状态用的是带符号的单一净功率值，符号正负代表放电或充电。这样做的好处是规避了同时充放电的约束。

约束校验我单独写了个函数，检查SOC递推是否越限、爬坡是否越限、P_grid是否越限。对越限的部分，我不是简单罚完就完事，而是优先做局部修复：比如把SOC越限的功率按比例缩回边界。修复后仍然不满意的，再在目标函数上叠加一个比较大的惩罚项。分层处理下来，种群里可行率提高非常多。

帕累托前沿提取和最优折中解选择，则在后处理阶段完成。具体实现是用非支配排序把所有个体分层，第一层即帕累托前沿；然后用模糊隶属度函数在所有非劣解里选折中解，也就是对每个解，分别计算它在两个目标上的归一化满意度再求和，取总满意度最高的那个作为最终的推荐调度方案。这个方法在实际论文和报告里接受度很高，因为避免了拍脑袋定权重。

4. 仿真结果分析与场景验证

4.1 帕累托前沿结果解读：经济与环保之间的取舍

跑完NSGA-II后，把所有非支配解画在二维坐标上，横轴是运行成本，纵轴是等效排放。我实际得到的帕累托前沿是一条单调递减的曲线，从左上到右下。左端点对应最低排放方案，成本最高；右端点对应最低成本方案，排放最高。中间解则分布在成本与排放的不同权衡上。

这个曲线的形状非常直观。成本最低的调度方案倾向于让柴油机在负荷高峰时段高负载运行，因为本地柴油机的综合发电成本可能低于峰时段从电网购电的价格，所以哪怕排放高，经济上也划算。排放最低的方案则反着来，尽可能压低柴油机出力，把光伏和风电尽量全部消纳，同时让储能在光伏大发时段充电、晚间放电，缺电时宁可高价购电也不让柴油机启动。

我选的折中解位于前沿中段，单位成本下降带来的排放上升幅度开始变大的那个拐点附近。这个位置对应的方案兼顾了两边：柴油机不频繁启停，储能有较大充放电深度，弃风弃光率控制在一个合理的范围内。

4.2 典型调度方案输出分析

把折中解对应的24小时出力曲线打出来看，典型场景是整个白天光伏出力旺盛，柴油机基本维持下限运行，储能持续充电；傍晚负荷开始爬升，光伏出力衰减，柴油机出力逐渐抬升；晚间负荷达到晚高峰，储能开始放电，配合柴油机一起顶峰；夜间风电出力较高，柴油机回落。

可以看出，储能的作用主要体现在“削峰填谷”：白天存光伏和风电的余量，晚间高峰释放，减少从电网买电的量。而柴油机则是兜底手段，既受爬坡约束限制，又在成本目标推动下尽可能优化使用。多台柴机之间还有一个负荷分配逻辑，我代码里采用了等微增率的思想做初始化偏好，不过因为NSGA-II的进化机制，它自己也能在种群中找到近似最优的分配方式。

此外我还加了两个验证场景：一个场景是柴油机和储能都正常运行，另一个场景是设置储能故障退出运行。结果很明显，储能退出后，系统更依赖柴油机和电网购电，运行成本上升，排放也上升，且晚间高峰时段的联络线功率经常触顶。这从侧面验证了模型约束的有效性。

5. 常见问题与排查技巧实录

5.1 高频问题速查表

我把实际编码过程中遇到的高频问题整理成一张表，方便大家直接对照排查：

5.2 参数整定与避坑经验

NSGA-II的参数看着不多，但整定起来还是有不少讲究。我常用的初值是：种群规模100，迭代200代，交叉概率0.9，变异概率取1/决策变量数，模拟二进制交叉的分布指数取20，多项式变异的分布指数取20。跑完看帕累托前沿的分布情况，如果前沿不完整，优先增加种群规模而不是迭代次数，因为迭代次数上去了，前沿可能还是稀疏，但种群规模越大，每个迭代代数上分布的个体也越多。

还有几个我在代码结构上的心得，常规博客里很少提到。

第一，目标函数里尽量少用全局变量。所有参数都塞进para结构体传进去，虽然会多打几个字，但调试时不会出现“某个变量在别的脚本里被改了”的诡异问题。

第二，决策变量边界别卡太死。比如柴油机最小出力取15kW，初始化时直接在15~50kW均匀采样，但爬坡约束只能靠下面几代的自然选择慢慢过滤，所以初代种群可行率不会太高。不用慌，把修复函数做扎实，几代之后可行率就上来了。

第三，MATLAB的数组运算比循环快得多，能用矩阵化写法就不要用几百次的for。尤其是计算两天内所有时段的燃料成本、SOC递推，用矩阵整体运算可以快一个数量级。SOC递推虽然看起来有自然的时序关系，实际上也可以写成向量化递推，只要写成SOC = SOC0 + cumsum(...)即可。

第四，如果想做进一步扩展，比如把确定性调度升级成考虑风光预测误差的随机优化，代码的扩展点主要集中在数据生成和目标函数两层。可以用蒙特卡洛生成多个风光场景，然后在目标函数里做期望值计算，算法部分基本不用动。

我在实际跑这个案例时，印象最深的是“残差变量法”带来的改变。一开始我是把所有变量都放进决策向量，结果跑了50代，有30代都在想办法平衡功率，帕累托前沿稀稀拉拉。改成电网功率作为残差反推后，同样的计算量，前沿质量和收敛速度都好了很多。这种“把等式约束藏进变量设计”的思路，其实在很多优化领域都通用，只是微网调度里的体现特别明显。

另外想提醒一句：网上MATLAB代码零零散散很多，但如果不能回答“为什么这里要这样解码”“为什么这个约束要用罚函数”，拿过来大概率也是跑不动的。真正花时间把模型和代码打通之后，你会觉得多目标优化调度其实没那么玄，它就是一个带约束的搜索问题，耐心调参，总能得到一张清晰漂亮的帕累托前沿。后续如果还有机会，我再写一版把风光预测不确定性和低碳运行约束加进去的扩展版本。