做微网优化的同行应该都有体会光搭一个纯电系统的优化模型练手半天就能跑通但一旦把燃气锅炉、储热罐、余热回收、CHP机组全部拉进来多能互补四个字马上就把问题复杂度抬高一个量级。热电联供型微网的核心矛盾在于电和热两套系统在物理上强耦合但在需求侧又各自独立变化冬天热负荷飙升时电负荷不一定大白天电负荷高峰时热负荷可能反而很低。这种电-热剪刀差正是优化调度要解决的主要矛盾。这篇文章我会完整拆解一套基于多能互补思路的热电联供型微网优化运行方案给出可复用的Matlab建模方法、约束方程写法和求解思路并附上我在实际调模型中踩过的坑。内容上不追求面面俱到但保证每一步都符合工程习惯适合正在做微网调度、综合能源系统优化、或者准备相关方向毕业设计的同学参考。1. 项目思路与优化问题怎么定1.1 热电联供微网里的多能互补到底指什么先捋清楚概念。我们说的热电联供型微网Combined Heat and Power Microgrid本质是一个同时承担电负荷与热负荷的小型能源系统常见组成包括燃气轮机或内燃机组成的CHP机组发电同时回收余热燃气锅炉作为备用热源电锅炉或热泵作为电转热环节实现以电补热储电、储热装置起到削峰填谷作用外部电网购电通道以及可能接入的光伏、风电等可再生能源多能互补体现在两个层面一是能源输入侧的多品种组合气、电、热、可再生能源协同供能二是转换环节的灵活性当热负荷紧缺时可以从电网购电驱动电锅炉补热当电价高、气价低时燃气机组多发电的同时多供热这就是互补的本质。很多新手在建模时容易把电系统和热系统拆开来独立做然后结果看似收敛实际不具备物理可行性。原因很简单CHP机组的电出力和热出力共用同一个原动机你不能让燃气轮机只发电不产热也不能只产热不发电。这种强耦合关系必须写进约束里否则优化结果就是空中楼阁。1.2 目标函数选最小运行成本还是最小碳排放优化运行的核心是目标函数。我见到的绝大多数研究都采用最小化系统日运行成本因为这是工程上最直接、最可解释的指标。典型成本项包括购电成本从上级电网购电的电量电费燃料成本燃气轮机和燃气锅炉的天然气消耗费用设备启停成本机组启停带来的寿命损耗费用运维成本按出力比例折算的设备维护费用也可以把碳排放量纳入做多目标优化。但我要提醒多目标并不是简单地把成本和排放加权碳定价法是工程上更合理的做法——给每吨CO2设定一个价格放进成本函数里统一优化。这样既避免了加权系数的主观性又能直接输出单一最优解。我常用的目标函数形式如下% 目标函数购电成本 燃料成本 碳成本 运维成本 objective sum(C_grid.*P_grid) ... % 购电成本 sum(C_gas.*F_chp C_gas.*F_boiler) ... % 燃料成本 sum(C_carbon.*(E_chp E_boiler)) ... % 碳排放成本 sum(K_om.*P_chp K_om.*P_grid); % 运维成本1.3 为什么选混合整数线性规划(MILP)框架优化运行问题的建模框架有很多种选择动态规划、遗传算法、粒子群、混合整数线性规划。我的建议是除非你有特殊需求否则优先选MILP。原因有三。第一MILP的全局最优性有数学保证智能算法只能给近似解在学术对比里容易被质疑第二MILP模型的约束可以直接对应物理逻辑调试直观第三现在的商业求解器对MILP的处理能力已经很强一台普通电脑上求解几百个变量的调度问题通常只需几十秒到几分钟。使用MILP的前提是目标函数和约束是线性的。这个条件在实际中基本可以满足成本曲线可以分段线性化设备效率可以取常数机组启停用0/1变量表示。光伏和负荷预测的波动可以通过场景法处理而场景法本身也只是把问题规模变大不改变MILP框架。2. 核心约束条件建模从物理规律到数学公式2.1 电功率平衡和热功率平衡两条线不能乱优化模型的最基本约束是功率平衡。简单说任何时刻系统发出的电功率必须等于电负荷加上损耗和储能充放电功率。写成方程就是% 电功率平衡约束T为调度时段数 Constraints [Constraints, P_grid P_chp P_pv P_dis P_load P_ch P_eb]; % 注意P_eb是电锅炉消耗的电功率热功率平衡类似所有热源产生的热功率必须覆盖热负荷% 热功率平衡约束 Constraints [Constraints, H_chp H_boiler H_dis H_load H_ch];这里有个细节容易出错储热装置也有充放热状态所以热平衡里同样要有储热充放变量。而且储电、储热不能同时充放这需要加约束% 储能不能同时充放M为足够大的数用于Big-M法 Constraints [Constraints, P_ch M.*u_ch, P_dis M.*u_dis, u_ch u_dis 1];2.2 CHP机组电热耦合可行域关键中的关键CHP机组是整个模型里最有技术含量的约束部分。燃气轮机的电出力和热回收之间不是简单的线性比例关系而是存在一个凸多边形可行域。对于典型的抽汽式CHP机组在电-热功率平面上可行运行区域可以描述为几个不等式围成的区域% CHP机组可行域约束以最小出力-最大出力-最大热出力围成的凸多边形为例 % 变量P_chp为电出力H_chp为热出力 Constraints [Constraints, P_chp max(P_min, 0.4*H_chp 15)]; % 最小电出力随热出力抬升 Constraints [Constraints, P_chp P_max - 0.2*H_chp]; % 最大电出力随热出力下降 Constraints [Constraints, H_chp H_max]; % 最大热出力限制这里的系数0.4和0.2来自机组的热电特性曲线实际工程数据需要通过机组厂家提供的工况图取值。可行域内的任意点代表一个可运行状态优化算法会在其中寻找最优工作点。要点是约束边界必须能反映发热时发电能力会受限这个物理规律。很多初稿模型只写了P_chp在一定区间内、H_chp在一定区间内但把两者写成相互独立的变量这会让优化结果出现又发电又大量发热的不合理工况。2.3 储能约束和机组爬坡约束别忘了时间耦合储能设备的约束要体现跨时段的状态转移。储电设备的电量SOC按如下方式递推for t 2:T Constraints [Constraints, SOC(t) SOC(t-1) P_ch(t)*eta_ch - P_dis(t)/eta_dis]; end % SOC要限制在允许范围 Constraints [Constraints, SOC_min SOC SOC_max]; % 初始和末尾SOC保持一致循环调度常用 Constraints [Constraints, SOC(1) SOC_init, SOC(T) SOC_init];爬坡约束刻画了机组出力调整的物理限制。燃气轮机每分钟能升多少负荷是有限的在15分钟粒度下有如下约束% 爬坡约束上调/下调速率限制 Constraints [Constraints, P_chp(t1) - P_chp(t) ramp_up]; Constraints [Constraints, P_chp(t) - P_chp(t1) ramp_down];我补充一个容易被忽略的细节如果做了机组启停优化还要加入启停相关的逻辑约束比如停机状态下出力为0和启动后最小运行时间限制。这一组约束才是MILP比连续LP复杂的地方。2.4 可再生能源出力的不确定性场景法处理很多论文会加入光伏和风电那么如何处理它们的间歇性就成了绕不开的问题。最常见做法是采用多场景随机优化生成若干个典型日的出力曲线给每个场景赋一个概率目标函数变成所有场景下的期望成本最小化。这样模型从确定性优化变成两阶段随机优化规模会成倍增长。如果不想一开始就碰随机优化我建议先跑确定性模型把光伏出力当成已知曲线来处理先跑通整个框架。等确定性版本稳定后再扩展为场景法。这个循序渐进的做法在硕士论文和实际工程项目里都适用可以避免一上来就被不确定性和求解规模困住。3. Matlab实现细节代码框架与核心模块3.1 建模工具怎么选YALMIP还是直接用求解器APIMatlab环境下做优化建模主流方案是YALMIP工具箱第三方求解器。YALMIP提供了高层的建模语法让你能像写数学表达式一样写约束和目标函数然后自动翻译成求解器可识别的格式。这样做的好处是代码结构清晰想换求解器只需改一行配置。当然也可以用linprog和intlinprog等原生函数但需要手动把约束写成矩阵形式。变量一多矩阵维度很容易错乱调试起来非常痛苦。除非你的模型极小否则我不推荐手写矩阵形式。求解器方面我实测过的几款如下表求解器许可证对MILP的处理能力适用场景Gurobi商业授权/学术免费很强大规模复杂问题首选CPLEX商业授权/学术免费很强传统老牌求解器COPT杉数商业授权/学术免费较强国产求解器文档友好SCIP开源免费中等小规模验证够用intlinprog内置Matlab自带够用小算例、快速演示我日常用Gurobi最多求解速度快数值稳定性好。如果是学生做课程设计先装SCIP免费版也能满足大部分场景。% 在Matlab中设置YALMIP调用Gurobi求解 options sdpsettings(solver,gurobi,verbose,1,showprogress,1);3.2 代码结构设计模块化拆分开别写成一坨我建议按下面的模块划分代码目录这是我从几个落地项目里沉淀下来的习惯data.m输入数据包括负荷曲线、能源价格、设备参数、预测曲线model_vars.m定义所有优化变量连续变量用sdpvar0/1变量用binvarconstraints.m按类别添加约束功率平衡、机组、储能、爬坡、备用objective.m构建目标函数run_optimization.m主程序依次调用上述模块并输出结果plot_results.m结果可视化变量命名保持类型_设备_时段的规则例如P_chp(t)、H_boiler(t)、SOC_ees(t)。这种命名方式在约束写多了之后极大降低排查成本。调度周期用1小时为步长一天24个时段起步如果做高精度调度可以细化到15分钟96个时段。每个时段的决策变量包括CHP电出力、CHP热出力、锅炉热出力、电锅炉电功率、储电充放、储热充放、购电功率等加起来大约10个变量乘以96个时段即约1000个变量其中一小部分是0/1变量。这个规模对现代求解器很轻松。3.3 一个最小的可运行模型直接套用下面给一个简化版模型核心代码。它以一天24小时为周期包含CHP机组、燃气锅炉、电锅炉、储能、购电通道目标是成本最小化。复制即用体验一下完整流程。%% 最小可运行的热电联供微网优化模型 T 24; P_load [数据导入]; % 电负荷曲线1x24 H_load [数据导入]; % 热负荷曲线1x24 Price [数据导入]; % 分时电价1x24 C_gas 3.4; % 天然气价格元/m3 % 1. 定义变量 P_chp sdpvar(1, T); % CHP电出力 H_chp sdpvar(1, T); % CHP热出力 P_boiler sdpvar(1, T); % 锅炉热出力kW P_eb sdpvar(1, T); % 电锅炉电功率 P_grid sdpvar(1, T); % 购电功率 SOC_e sdpvar(1, T); % 储电SOC SOC_h sdpvar(1, T); % 储热SOC u_chp binvar(1, T); % CHP启停变量 % 2. 约束集合 Constraints []; % 电功率平衡 for t 1:T Constraints [Constraints, P_grid(t) P_chp(t) P_load(t) P_eb(t)]; end % 热功率平衡CHP热出力 锅炉补热 热负荷 电锅炉耗热? 实际电锅炉耗电生热 % 注意电锅炉产热功率记入热平衡本简化模型中H_eb COP * P_eb for t 1:T Constraints [Constraints, H_chp(t) P_boiler(t) COP_eb*P_eb(t) H_load(t)]; end % CHP可行域 for t 1:T Constraints [Constraints, P_chp(t) 20.*u_chp(t) 0.4.*H_chp(t)]; Constraints [Constraints, P_chp(t) 100.*u_chp(t) - 0.2.*H_chp(t)]; Constraints [Constraints, H_chp(t) 50.*u_chp(t)]; end % 储能SOC递推和容量限制 for t 2:T Constraints [Constraints, SOC_e(t) SOC_e(t-1) 0.9*P_ch_e(t) - P_dis_e(t)/0.9]; Constraints [Constraints, SOC_h(t) SOC_h(t-1) 0.95*P_ch_h(t) - P_dis_h(t)/0.95]; end Constraints [Constraints, 0.2 SOC_e 0.9, 0.2 SOC_h 0.9]; % 电网购电上限 Constraints [Constraints, 0 P_grid 200]; % 3. 目标购电成本 燃料成本CHP和锅炉 F_chp (P_chp 0.5*H_chp) ./ 9.7; % 天然气消耗量m3/h F_boiler P_boiler ./ (0.9 * 9.7); % 锅炉气耗 objective sum(Price.*P_grid) C_gas*(sum(F_chp) sum(F_boiler)); % 4. 求解 options sdpsettings(solver,gurobi,verbose,1); optimize(Constraints, objective, options);这个模型省略了储能充放电变量的定义但约束框架是完整可扩展的。想完整跑通的话把储能的充放变量补上再补齐SOC递推逻辑即可。这是最小可行的模板扩展性好。3.4 结果可视化学会读优化结果运行完模型后第一件事是画负荷平衡图。优化结果必须满足每个时段的功率平衡而画图能直观看出电负荷由谁供给、热负荷由谁承担。我常用两张图一张堆叠面积图展示电功率构成另一张展示热功率构成。figure; area(1:T, [P_grid, P_chp, P_eb], LineWidth, 1); legend(购电,CHP发电,电锅炉耗电,Location,best); xlabel(时段/h); ylabel(功率/kW); title(电功率平衡构成); saveas(gcf, electric_balance.png);热功率图同理。第三张图建议画SOC曲线看储电储热是否处在合理区间。好的调度结果SOC曲线应该是平滑的锯齿形峰谷时段有充放动作而不是一条直线或者剧烈振荡。4. 常见问题与排查技巧实录4.1 求解器加载失败或许可证报错YALMIP调用外部求解器时最常见的报错是“No suitable solver installed”或者是License错误。我的建议装求解器的时候务必把许可证路径配置好。Gurobi在Windows下通常安装到C盘YALMIP会自动识别但Mac和Linux下可能需要在MATLAB里手动执行% 手动添加Gurobi路径 addpath(C:\gurobi1100\matlab); % 以实际安装路径为准 gurobi_setup() % 初始化尽调如果是intlinprog自带的求解器则是许可证随Matlab没这个烦恼。但如果模型变量较多自带求解器速度会明显下滑。4.2 模型不可行怎么定位是哪个约束引发的不可行是MILP建模最容易踩的坑。解决办法是逐个约束加松弛变量% 把松弛量加到可能出问题的约束上 s sdpvar(1, T); % 非负松弛变量 Constraints [Constraints, P_grid(t) P_chp(t) s(t) P_load(t) P_eb(t)]; % 然后检查最优解里的s值哪个时段的s不为0问题就出在哪个时段把目标函数暂时改为minimize sum(s)求解后观察哪一时段松弛量最大就可以反过来检查该时段是负荷高得离谱还是机组出力上限设太小。我用这个方法排查过很多次模型bug比盲猜高效得多。4.3 几个很值得留意的实操细节第一所有变量的单位务必统一。我做过一个项目CHP发热量的单位用了GJ而电出力的单位用了kWh结果热平衡约束怎么检查都不收敛。统一成kW功率作为基准单位每条约束都代入单位复核一遍能省去大量排查时间。第二Big-M法里的M不能取太大。虽然理论上M取无穷大就行但数值上M过大会造成病态矩阵让求解器精度崩塌。M只取该变量物理上限的10倍左右数值上已经足够。第三供热网络本身的热惯性在标准模型里一般忽略。如果你想做更精细的控制则需要考虑热网管道蓄热效应这时不能再用静态热平衡要引入动态热网方程。第四目标函数里的平方项能用线性逼近就尽量线性。如果不小心写了二次项模型会变成MIQP求解难度立刻上一个台阶。除非确实需要否则保持线性。4.4 关于数据和场景的两点经验补充运行调度模型数据质量比模型本身更影响结果。我建议从公开数据平台抓取典型日数据比如美国能源信息署(EIA)的部分公开负荷数据或者国内高校开源数据仓库。自行构造假数据做演示可以但如果后续要对标真实场景务必找有据可查的数据源。对于场景法处理不确定性我建议先用“聚类法”把365天的数据聚成3-5个典型场景每个场景对应一个概率。这样既保留了主要天气类型又不至于让问题规模失控。随机优化在论文里的说服力比确定性优化明显更高而且聚类实现起来不难。5. 完整扩展方向的几条参考思路前面把主干流程都讲完了如果你还想把这套模型继续延展我给出三个比较现实的方向。第一是加入需求响应机制。通过引入可转移负荷和可削减负荷让电负荷曲线从刚性的变成柔性的。做法是在约束里增加负荷转移变量和补偿成本模型规模多出几十个变量但能明显提升微网运行的经济性。第二是引入热网管道动态模型。传统做法只做能源站内的调度如果把管网的传输延迟、热损失、储热特性也建模系统灵活性会再上一个台阶。这对大区域供热系统尤其有效代价是模型的微分方程更多求解时间变长。第三是探索日前-日内两阶段调度。日前用场景法做随机优化确定机组启停和储能计划日内用滚动优化跟踪实际负荷偏差。这种“先计划后调整”的结构是工程落地的主流做法也比单阶段模型更有说服力。这几个方向上我踩过不少坑比如日前计划用太激进会导致日内滚动优化频繁反推、热负荷预测不准导致储热罐放空等。实际工程里计划一定要留裕度运行边界宁可窄一点也不要在临界点跳舞。6. 写给自己也算写给正在做这个方向的你我从第一次接触微网优化到现在最深的感受是优化模型的建成只是第一步真正的功夫在“约束写得好不好”“结果能不能落回物理现场”。你写出的每一条不等式都应该能对应到现场设备的某个调节阀或者某个功率上限否则这个模型就只是个数学游戏。热电联供微网优化运行这个课题Matlab代码实现其实只是一个载体核心还是对能源转换物理过程的理解深度。建议你做任何改动时都反问自己如果把这个约束删掉现场会不会出问题如果把这个参数放大一倍设备会不会扛不住带着这种工程直觉去写模型才能做出真正有价值的东西。最后分享一个小技巧每次跑完优化把结果里的调度曲线导出来和实际运行曲线叠在一起看不仅是看趋势像不像还要看每个设备在峰值时有没有超出预期出力。有一次我发现CHP机组在优化结果中频繁启停而实际根本不允许这么操作就是加了最小开关机时间约束才解决问题。这类问题建模时想不起来等对图横竖不对的时候自然而然就想到了。