简介：本资源面向能源系统优化、综合能源调度方向的学习者与研究人员，提供一套基于微型燃气轮机的冷热电联供系统优化运行程序。程序选取夏季典型日，仅考虑电负荷与冷负荷，以燃气轮机出力和空调出力为决策变量，余热驱动吸收式制冷，电制冷协同参与，电量由燃机提供并可向电网购电，目标函数为联供系统全天运行费用最小，采用粒子群算法求解最优调度方案。压缩包共8个文件，约61KB，以m脚本文件为主，辅以运行结果图片和asv备份文件，代码每行均有注释，结构简练，便于理解算法实现与建模思路。目前已有253人学习下载。读者可获得完整的PSO求解流程、目标函数与约束建模方法、典型日调度结果及功率曲线图，适合作为课程设计、论文复现或算法入门的参考素材。
1. 微型燃气轮机冷热电联供优化：一份能跑出收敛曲线的 PSO 实战代码
夏季典型日的园区微网里，电负荷和冷负荷同时压上来，燃气轮机出力、吸收式制冷、电制冷、电网购电四条路径怎么分配，直接决定一天运行费用是三千还是五千。这份资源就是冲着这个场景来的：用粒子群算法求解微型燃气轮机冷热电联供系统的一天最优调度，目标函数是运行费用最小，约束里考虑了燃机出力上限、空调出力、余热回收给吸收式制冷供热，电量不足部分从电网购入。代码包里有initpop.m、calobjvalue.m、calvalue.m、opsomian.m这几个核心文件，外加冷功率曲线图和电功率曲线图两张结果图。适合正在做综合能源系统调度、需要一套能直接跑通并改参数的 PSO 代码的从业者，也适合拿它当教学案例拆解粒子群在机组组合问题上的落地方式。每一行都有注释，代码简练，出图效果直观，收敛速度在同类问题里属于能接受的水平。
2. 冷热电联供与 PSO 的匹配逻辑：为什么选它、模型怎么建
2.1 三联供系统的能量流与决策变量
冷热电联供的核心是把燃气轮机发电后的高温余热回收，驱动吸收式制冷机产冷，不足的冷负荷再由电制冷补足，电量缺口从电网买。这个系统的决策变量通常包括：燃气轮机逐时出力、吸收式制冷机逐时出力、电制冷机逐时出力、电网购电功率。夏季典型日只有电负荷和冷负荷，没有热负荷，所以余热全部导向吸收式制冷，模型比冬季三负荷场景简单，但约束耦合依然存在——燃机出力决定了余热上限，余热上限决定了吸收式制冷上限，冷负荷缺口再决定电制冷和购电。
常见做法是把一天分成 24 个时段，每个时段一组决策变量，粒子维度就是 24 乘以变量个数。目标函数写成一整天燃料成本加购电成本，燃料成本由燃机出力除以发电效率再乘以气价得到，购电成本由购电功率乘以分时电价得到。约束用罚函数处理，越限就加一个大惩罚项，这样 PSO 在迭代中会自动避开不可行域。
2.2 粒子群算法在机组调度上的适配性
PSO 处理这类连续变量、带约束、非线性目标的问题有天然优势：不需要梯度信息，种群并行搜索，收敛速度比遗传算法快，参数少，调起来直观。标准 PSO 的速度更新公式里，惯性权重控制全局探索与局部开发的平衡，学习因子 c1 和 c2 分别拉向个体最优和全局最优。在这份代码里，opsomian.m是主循环，initpop.m负责初始化种群，calobjvalue.m算目标函数值，calvalue.m做约束检查和罚函数计算。
选 PSO 而不是动态规划或混合整数规划，理由很实际：动态规划在 24 时段、多变量、连续出力场景下状态空间爆炸；混合整数规划需要线性化，建模工作量大，而且燃机效率曲线往往是非线性的。PSO 直接搜连续空间，代码量小，改目标函数和约束都方便，适合快速验证调度策略。代价是它不保证全局最优，但工程上一天费用差几个百分点通常可以接受，收敛快才是刚需。
2.3 目标函数与约束的代码落地
打开calobjvalue.m，核心逻辑是逐时段累加费用。下面这段是典型写法，我按资源里的结构还原并加了注释：
function [objvalue] = calobjvalue(pop, load_e, load_c, price_gas, price_e) % pop: 种群矩阵，每行一个粒子，列按 [燃机出力, 吸收式制冷, 电制冷, 购电] 逐时排列 % load_e: 24小时电负荷，load_c: 24小时冷负荷 % price_gas: 气价，price_e: 24小时分时电价 [N, ~] = size(pop); objvalue = zeros(N, 1); for i = 1:N x = pop(i, :); % 拆出各变量，每24个一组 P_gt = x(1:24); % 燃气轮机逐时出力 P_ac = x(25:48); % 吸收式制冷逐时出力 P_ec = x(49:72); % 电制冷逐时出力 P_buy = x(73:96); % 电网购电 % 燃料成本：燃机出力 / 发电效率 * 气价，效率取0.3 cost_fuel = sum(P_gt / 0.3 .* price_gas); % 购电成本 cost_buy = sum(P_buy .* price_e); objvalue(i) = cost_fuel + cost_buy; end逻辑说明：每个粒子被拉成 96 维向量，前 24 维是燃机出力，接着 24 维是吸收式制冷，再 24 维是电制冷，最后 24 维是购电。燃料成本用燃机出力除以发电效率再乘气价，购电成本直接乘分时电价。参数方面，发电效率 0.3 是微型燃机的典型值，气价和电价从外部传入，方便改。如果要把余热回收效率也建模进去，可以在P_ac的计算里加一个余热利用系数，但资源里为了简练没有展开。
约束检查在calvalue.m里，主要看三条：燃机出力不超过额定上限，吸收式制冷不超过余热可提供的冷量上限，电制冷加吸收式制冷等于冷负荷，燃机加购电等于电负荷。违反任一条就加罚。罚函数系数一般取目标函数量级的 10 到 100 倍，太小约束压不住，太大种群多样性丢失快。
2.4 初始化与主循环的参数设置
initpop.m负责生成初始种群。常见做法是在变量上下限之间均匀随机，但为了加速收敛，可以按负荷比例给一个初始猜测，再叠加扰动。资源里的初始化应该是均匀随机，因为代码简练优先。
opsomian.m是主循环，关键参数有：种群规模 N、最大迭代次数 maxgen、惯性权重 w、学习因子 c1 和 c2。下面是一个可抄的骨架：
% PSO主循环参数 N = 50; % 种群规模，太小易早熟，太大费时间 maxgen = 200; % 迭代次数，看收敛曲线定 w = 0.8; % 惯性权重，0.4~0.9之间调 c1 = 1.5; % 个体学习因子 c2 = 1.5; % 全局学习因子 % 初始化 pop = initpop(N, 96, lb, ub); v = zeros(N, 96); pbest = pop; pbest_obj = calobjvalue(pop, load_e, load_c, price_gas, price_e); [gbest_obj, idx] = min(pbest_obj); gbest = pop(idx, :); for gen = 1:maxgen for i = 1:N v(i,:) = w*v(i,:) + c1*rand*(pbest(i,:)-pop(i,:)) + c2*rand*(gbest-pop(i,:)); pop(i,:) = pop(i,:) + v(i,:); % 边界处理 pop(i,:) = max(pop(i,:), lb); pop(i,:) = min(pop(i,:), ub); % 计算目标并更新个体最优 obj = calobjvalue(pop(i,:), load_e, load_c, price_gas, price_e); if obj < pbest_obj(i) pbest(i,:) = pop(i,:); pbest_obj(i) = obj; end end [min_obj, idx] = min(pbest_obj); if min_obj < gbest_obj gbest_obj = min_obj; gbest = pbest(idx, :); end % 记录收敛曲线 trace(gen) = gbest_obj; end参数说明：种群规模 50 是中小规模问题的常用值，96 维变量下 50 个粒子能覆盖得比较开。迭代 200 次通常能看到收敛，如果曲线还在降就加到 500。惯性权重 0.8 偏探索，后期可以线性降到 0.4 增强开发。c1 和 c2 取 1.5 是经典配置，也可以试 2.0 和 2.0。边界处理用截断法，简单有效，缺点是粒子撞墙后速度没清零，可能来回震荡，改进做法是撞墙时速度反向或置零。
3. 从零跑通代码：环境、数据与出图
3.1 文件结构与运行顺序
资源包解压后能看到这些文件：initpop.asv、calvalue.m、initpop.m、冷功率曲线图.jpg、calobjvalue.m、opsomian.m、电功率曲线图.jpg、opso.asv。.asv是 MATLAB 的自动保存文件，不是主程序，真正要跑的是opsomian.m。运行顺序是：先确保initpop.m、calobjvalue.m、calvalue.m在同一目录，然后在 MATLAB 里直接运行opsomian.m。如果报未定义函数，检查当前工作目录是不是设到了代码所在文件夹。
冷功率曲线图和电功率曲线图是运行结果，横轴是 24 小时，纵轴是功率。冷功率图里能看到吸收式制冷和电制冷的分配，电功率图里能看到燃机出力和购电的分配。这两张图是验证调度方案是否合理的直接依据。
3.2 负荷与价格数据的准备
代码里负荷和价格数据通常是硬编码或从外部文件读。如果资源里是硬编码，找到opsomian.m开头的数据段，按下面格式替换成自己的数据：
% 24小时电负荷(kW) load_e = [120 115 110 105 110 130 160 190 210 220 225 220 ... 215 220 225 230 235 240 230 210 190 170 150 130]; % 24小时冷负荷(kW) load_c = [80 75 70 70 75 90 120 150 180 200 210 205 ... 200 205 210 215 220 225 215 195 170 150 120 95]; % 分时电价(元/kWh)，峰平谷三段 price_e = [0.4 0.4 0.4 0.4 0.4 0.4 0.8 0.8 1.2 1.2 1.2 1.2 ... 1.2 1.2 1.2 1.2 0.8 0.8 0.8 0.8 0.4 0.4 0.4 0.4]; % 气价(元/kWh) price_gas = 0.35;逻辑说明：电负荷和冷负荷按夏季典型日双峰曲线给，早晚高峰明显。电价用峰平谷三段，峰段 1.2 元，平段 0.8 元，谷段 0.4 元，这样 PSO 会自动倾向于在谷段多购电、峰段少购电。气价 0.35 元是常见工业气价量级。参数改完后直接跑，收敛曲线和两张功率图会更新。
3.3 收敛曲线与结果图的判读
跑完后 MATLAB 会弹出收敛曲线，横轴迭代次数，纵轴最优费用。正常形态是前 50 代快速下降，100 代后趋于平缓，200 代基本不变。如果曲线一直震荡不降，检查罚函数系数是不是太大导致可行域被压扁，或者种群规模太小。如果曲线在 20 代就平了，可能是早熟，把惯性权重调大或种群规模加到 80。
冷功率曲线图里，吸收式制冷应该跟随余热可用量走，电制冷补缺口。如果电制冷全天顶满而吸收式很低，说明燃机出力被压得太低，可能是气价设太高或电价设太低。电功率图里，购电应该在峰段明显下降，谷段上升，这是分时电价在起作用。如果购电曲线平的，检查电价数组是不是没生效。
3.4 改目标函数与加约束的实操
想把目标函数从费用最小改成费用加碳排放最小，在calobjvalue.m里加一项：
% 碳排放成本：购电碳排放 + 燃机碳排放 carbon_buy = sum(P_buy * 0.8); % 购电碳排放因子0.8 kg/kWh carbon_gas = sum(P_gt / 0.3 * 0.2); % 燃机碳排放因子0.2 kg/kWh cost_carbon = (carbon_buy + carbon_gas) * 0.1; % 碳价0.1元/kg objvalue(i) = cost_fuel + cost_buy + cost_carbon;参数说明：购电碳排放因子 0.8 是煤电为主的电网典型值，燃机碳排放因子 0.2 是天然气发电的常见值，碳价 0.1 元每千克是示例，实际按当地碳市场价改。加约束的话，在calvalue.m里增加判断，比如燃机爬坡率限制，相邻时段出力差不超过额定出力的 20%，违反就加罚。
4. 避坑与排查：跑不通、不收敛、结果反直觉的常见原因
4.1 现象：运行报错“未定义函数或变量”
原因：MATLAB 当前工作目录不在代码所在文件夹，或者.m文件名和函数名不一致。资源里opsomian.m是主脚本，initpop.m、calobjvalue.m、calvalue.m是函数文件，函数名必须和文件名一致。
解决：在 MATLAB 里用cd切到代码目录，或者把目录加到路径。检查每个函数文件第一行function后面的名字和文件名是否完全一致，大小写敏感。
4.2 现象：收敛曲线剧烈震荡，费用不降
原因：罚函数系数过大，可行域被压成一条线，粒子一进去就被弹飞；或者惯性权重太小，粒子飞不动。
解决：先把罚函数系数降到目标函数量级的 10 倍左右试，再逐步加。惯性权重从 0.8 开始，如果震荡就降到 0.6，如果收敛太慢就升到 0.9。学习因子 c1 和 c2 不要同时大于 2.0，否则速度爆炸。
4.3 现象：结果图里购电全天为零或全天顶满
原因：电价数组没传进去，或者购电上下限设成了固定值。也有可能是燃机出力上限设得太大，PSO 发现全用燃机更便宜，购电自然为零。
解决：检查price_e是不是 24 维，检查购电变量上下限是不是[0, 额定购电上限]。如果燃机气价设得太低，确实会全用燃机，这时候看冷负荷能不能被余热满足，满足不了就必须购电或电制冷，检查约束有没有生效。
4.4 现象：吸收式制冷出力始终为零
原因：余热计算没接上，或者吸收式制冷的上下限设成了零。资源里余热来自燃机，如果燃机出力被压得很低，余热不够，吸收式自然为零。
解决：检查calvalue.m里余热和吸收式制冷的关联约束，确认P_ac的上限是余热乘以制冷系数。如果燃机出力本身很低，看目标函数里气价是不是太高，把气价降到 0.2 试一下，燃机出力会上去，余热和吸收式跟着起来。
4.5 现象：每次运行结果不一样，费用差很多
原因：PSO 是随机算法，初始种群和速度更新都带随机数，没有固定随机种子。这是正常现象，但差太多说明算法不稳定。
解决：在opsomian.m开头加rng(1)固定随机种子，方便复现。如果固定种子后不同种子结果差异大，说明种群规模不够，加到 80 或 100。也可以跑多次取最优，工程上常见做法是跑 10 次留最好的那组调度方案。
5. 进阶技巧：把 PSO 调稳、把调度方案用起来
5.1 惯性权重线性递减与收敛判据
标准 PSO 的惯性权重固定，前期探索后期开发的需求它兼顾不了。我一般会改成线性递减：
w_start = 0.9; w_end = 0.4; for gen = 1:maxgen w = w_start - (w_start - w_end) * gen / maxgen; % 后续速度更新用这个w end这样前 50 代 w 在 0.9 附近，粒子敢飞，全局搜索强；后 100 代 w 降到 0.5 以下，粒子收着飞，局部开发细。收敛判据不要只看迭代次数，加一个早停：连续 30 代最优值变化小于 0.1%，就跳出，省时间。这个技巧在 96 维问题上能省 30% 到 50% 的迭代。
5.2 约束处理从罚函数换成可行性规则
罚函数调参麻烦，系数敏感。另一种常见做法是可行性规则：比较两个粒子时，可行粒子永远优于不可行粒子；都可行比目标函数，都不可行比违反约束的程度。这样不用设罚函数系数，但代码要改比较逻辑。在opsomian.m里更新个体最优和全局最优的地方，把直接比目标值改成先比可行性标志。资源里用的是罚函数，适合入门，进阶可以试可行性规则。
5.3 把调度方案导出成表格给运行人员
PSO 跑完得到gbest向量，把它拆成 24 小时的燃机出力、吸收式制冷、电制冷、购电，写进 Excel 或 CSV：
P_gt = gbest(1:24); P_ac = gbest(25:48); P_ec = gbest(49:72); P_buy = gbest(73:96); T = table((1:24)', P_gt', P_ac', P_ec', P_buy', ... 'VariableNames', {'Hour','GasTurbine','AbsorptionChiller','ElectricChiller','GridPurchase'}); writetable(T, 'dispatch_plan.csv');运行人员拿到这张表，每小时该开多少燃机、多少吸收式、买多少电，一目了然。注意单位统一成 kW，如果原始数据是 MW，导出前乘 1000。这张表也是和实际运行数据对比的依据，跑一段时间后把实际值和计划值画在一起，看偏差在哪几个时段大，反过来调负荷预测或模型参数。
5.4 多目标扩展与场景切换
夏季典型日只有电和冷，冬季要加熱负荷，过渡季可能只有电。切换场景时，改load_c和load_e的数组，同时在calvalue.m里调整约束：冬季余热优先供热，吸收式制冷可以关掉或改成供热。多目标的话，把费用和碳排放两个目标加权成单目标，权重从 0.1 到 0.9 扫一遍，得到帕累托前沿的近似，挑一个折中点。我一般会跑 9 个权重组合，每个跑 3 次取最优，画一条费用-碳排放曲线，给决策者看。
从那以后我每次拿到新的综合能源调度代码，都强制先跑一遍原始数据，把收敛曲线和两张功率图存下来当基线，再改参数。改完对比基线，如果费用降了但某台设备出力全天顶满，就回去查约束是不是漏了。这套流程帮我省了不少返工。希望帮到你。
本文还有配套的精品资源，点击获取
